How to solve a boundary value problem x"+x'=1, x(0)=1, x(1)=2 , with an impulse condition x(i/3+)=2x(i/3)+1, x'(i/3+)=3x'(i/3)-1;

Question

Nauryzbay 2023-6-22

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1986814-how-to-solve-a-boundary-value-problem-x-x-1-x-0-1-x-1-2-with-an-impulse-condition-x-i-3-2

编辑： Torsten 2023-6-23

function impulsbound %x"+x'=1, x(0)=1, x(i/3+)=2*x(i/3)+1, x'(i/3+)=3*x'(i/3)-1, i=1,2. x(1)=2;

format long

clear all

t=linspace(0,1,100);

solinit = bvpinit(t,[1 2]) ;

x1=1;

u=2;

teta(1)=0;

%teta(u+2)=0.2;

for i=1:u

teta(i+1)=i/3;

end

for j=1:1

for k=1:u

%initial_func=[x1,x2];

%[t,x] = ode45(@s,[teta(k):0.0001:teta(k+1)],[x1(j) x2(j)]);

x = bvp4c(@s,@bc2,solinit,[teta(k):0.0001:teta(k+1)],[x1(j) x2(j)]);%

n=length(t);

%disp(size(x));

x1(j)=x(n,1)+x(n,1)+1;

x2(j)=x(n,2)+2*x(n,2)-1;

hold on

%view(30,15);

hold on

figure(1)

subplot(2,1,1);

plot(t,x(:,1),'color','r','Linewidth',1.2);

xlabel('\bf t'); ylabel('$$z$$','interpreter','latex','fontsize',16); zlabel('\bf \psi_2');

grid on

hold on

subplot(2,1,2);

plot(t,x(:,2),'color','r','Linewidth',1.2);

xlabel('\bf t'); ylabel('$$y$$','interpreter','latex','fontsize',16); zlabel('\bf \psi_3');

grid on

hold on

end

function dx=s(t,x)

dx=zeros(2,1); % создает нулевой вектор-столбец

dx(1)=x(2);

dx(2)=1-x(1);

function res = bc2(xa, xb)

res = [xa(1)-1

xb(1)-2];

Reference to a cleared variable x2.

Error in impulsbound (line 17)

x = bvp4c(@s,@bc2,solinit,[teta(k):0.0001:teta(k+1)],[x1(j) x2(j)]);%

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Torsten 2023-6-22

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1986814-how-to-solve-a-boundary-value-problem-x-x-1-x-0-1-x-1-2-with-an-impulse-condition-x-i-3-2#answer_1260829

编辑：Torsten 2023-6-22

在 MATLAB Online 中打开

syms t x1(t) x2(t)

eqn1 = diff(x1,t,2) + diff(x1,t) == 1;

Dx1 = diff(x1,t);

eqn2 = diff(x2,t,2) + diff(x2,t) == 1;

Dx2 = diff(x2,t);

eqns = [eqn1,eqn2];

conds = [x1(0) == 1, x2(1/3) == 2*x1(1/3)+1, Dx2(1/3) == 3*Dx1(1/3)-1, x2(1) == 2];

sol = dsolve(eqns,conds)

sol = struct with fields:

x2: t + (2*exp(1) + 2*exp(1/3) + exp(2/3) - 12)/(2*exp(1) + exp(2/3) - 3) - (exp(-t)*exp(1/3)*(2*exp(1) - 9*exp(2/3)))/(2*exp(1) + exp(2/3) - 3) x1: t - (4*exp(1) - 3*exp(1/3) - 3*exp(2/3) + 9)/(3*(2*exp(1) + exp(2/3) - 3)) + (exp(-t)*(10*exp(1) - 3*exp(1/3)))/(3*(2*exp(1) + exp(2/3) - 3))

hold on

fplot(sol.x1,[0,1/3])

fplot(sol.x2,[1/3,1])

hold off

grid on

3 个评论
显示 1更早的评论隐藏 1更早的评论

Torsten 2023-6-22

编辑：Torsten 2023-6-23

在 MATLAB Online 中打开

Then you will have to write a little more: you will have equations and transmission conditions for x1(t),...,x10(t).

And also the above code is set up for only one transmission condition at x=1/3, not both at x=1/3 and x=2/3. So already in this case, you have to use x1(t), x2(t) and x3(t) as follows:

syms t x1(t) x2(t) x3(t)

eqn1 = diff(x1,t,2) + diff(x1,t) == 1;

Dx1 = diff(x1,t);

eqn2 = diff(x2,t,2) + diff(x2,t) == 1;

Dx2 = diff(x2,t);

eqn3 = diff(x3,t,2) + diff(x3,t) == 1;

Dx3 = diff(x3,t);

eqns = [eqn1,eqn2,eqn3];

conds = [x1(0) == 1, ...

x2(1/3) == 2*x1(1/3)+1, Dx2(1/3) == 3*Dx1(1/3)-1,...

x3(2/3) == 2*x2(2/3)+1, Dx3(2/3) == 3*Dx2(2/3)-1,...

x3(1) == 2];

sol = dsolve(eqns,conds)

sol = struct with fields:

x2: t - (10*exp(1) - 57*exp(1/3) + 5*exp(2/3) + 108)/(3*(4*exp(1) + 3*exp(1/3) + 2*exp(2/3) - 9)) - (exp(-t)*exp(1/3)*(4*exp(1) + 3*exp(1/3) - 29*exp(2/3)))/(4*exp(1) + 3*exp(1/3) + 2*exp(2/3)… x1: t - (19*exp(1) - 18*exp(1/3) - 6*exp(2/3) + 27)/(3*(4*exp(1) + 3*exp(1/3) + 2*exp(2/3) - 9)) + (exp(-t)*(31*exp(1) - 9*exp(1/3)))/(3*(4*exp(1) + 3*exp(1/3) + 2*exp(2/3) - 9)) x3: t + (4*exp(1) + 17*exp(1/3) + 6*exp(2/3) - 93)/(4*exp(1) + 3*exp(1/3) + 2*exp(2/3) - 9) - (2*exp(-t)*exp(2/3)*(2*exp(1) - 42*exp(1/3) + 7*exp(2/3)))/(4*exp(1) + 3*exp(1/3) + 2*exp(2/3) - 9)

hold on

fplot(sol.x1,[0,1/3])

fplot(sol.x2,[1/3,2/3])

fplot(sol.x3,[2/3,1])

hold off

grid on