How to find the jacobian of the system of nonlinear equations?

1 Ansicht (letzte 30 Tage)
I wrote a following program, but it is giving some errors?
clear all
syms x
h=4;
gma=18.4;
ka1=0.2;
kp1=8.76
ka2=0.2;
kp2=8.76;
sindel=0.437;
cosdel=0.9;
pa1=ka1*gma*(h*x(1)+0.5*x(1)^2);
pp1=kp1*gma*0.5*x(1)^2;
pa2=ka2*gma*(x(1)*x(2)+0.5*x(2)^2);
pp2=kp2*(x(2)*(h+x(1))+0.5*x(2)^2);
za1=(0.5*h*x(1)^2+(x(1)^3)/6)/(h*x(1)+0.5*x(1)^2);
zp2=(0.5*(h+x(1))*x(2)^2+((x(2)^3)/3))/((h+x(1))*x(2)+0.5*x(2)^2);
za2=(0.5*x(1)*x(2)^2+((x(2)^3)/3))+(x(1)*x(2)+0.5*x(2)^2);
e1=pp1*sindel-pa1*sindel-pp2*sindel-pa2*sindel;
e2=pp1*cosdel+pa2*cosdel-pa1*cosdel-pp2*cosdel-x(3);
e3=pp1*cosdel*(x(1)/3)+pp2*cosdel*zp2-pa1*cosdel*za1-pa2*cosdel*za2-x(3)*(x(1)+(h/3));
g=[e1; e2; e3];
J=jacobian([e1, e2, e3], [x(1), x(2), x(3)])

Akzeptierte Antwort

Stephan
Stephan am 23 Jan. 2019
Hi,.
try:
syms x1 x2 x3
h=4;
gma=18.4;
ka1=0.2;
kp1=8.76
ka2=0.2;
kp2=8.76;
sindel=0.437;
cosdel=0.9;
pa1=ka1*gma*(h*x1+0.5*x1^2);
pp1=kp1*gma*0.5*x1^2;
pa2=ka2*gma*(x1*x2*+0.5*x2^2);
pp2=kp2*(x2*(h+x1)+0.5*x2^2);
za1=(0.5*h*x1^2+(x1^3)/6)/(h*x1+0.5*x1^2);
zp2=(0.5*(h+x1)*x2^2+((x2^3)/3))/((h+x1)*x2+0.5*x2^2);
za2=(0.5*x1*x2^2+((x2^3)/3))+(x1*x2+0.5*x2^2);
e1=pp1*sindel-pa1*sindel-pp2*sindel-pa2*sindel;
e2=pp1*cosdel+pa2*cosdel-pa1*cosdel-pp2*cosdel-x3;
e3=pp1*cosdel*(x1/3)+pp2*cosdel*zp2-pa1*cosdel*za1-pa2*cosdel*za2-x3*(x1+(h/3));
g=[e1; e2; e3];
J=jacobian([e1, e2, e3], [x1, x2, x3])
results in:
J =
[ (1075457*x1)/15625 - (95703*x2)/25000 - (10051*x2^3)/12500 - 20102/3125, - (95703*x1)/25000 - (95703*x2)/25000 - (30153*x1*x2^2)/12500 - 95703/6250, 0]
[ (207*x2^3)/125 - (1971*x2)/250 + (88596*x1)/625 - 1656/125, (621*x1*x2^2)/125 - (1971*x2)/250 - (1971*x1)/250 - 3942/125, -1]
[ (44298*x1^2)/625 - x3 - (207*x2^3*(x1*x2 + (x1*x2^2)/2 + x2^2/2 + x2^3/3))/125 - (1656*x1)/125 + (x2^2*((1971*x2*(x1 + 4))/250 + (1971*x2^2)/500))/(2*(x2^2/2 + (x1 + 4)*x2)) - (((414*x1)/125 + 1656/125)*(x1^3/6 + 2*x1^2))/(x1^2/2 + 4*x1) + (1971*x2*(x2^3/3 + (x1/2 + 2)*x2^2))/(250*(x2^2/2 + (x1 + 4)*x2)) - (207*x1*x2^3*(x2^2/2 + x2))/125 + (((207*x1^2)/125 + (1656*x1)/125)*(x1^3/6 + 2*x1^2)*(x1 + 4))/(x1^2/2 + 4*x1)^2 - (x2*(x2^3/3 + (x1/2 + 2)*x2^2)*((1971*x2*(x1 + 4))/250 + (1971*x2^2)/500))/(x2^2/2 + (x1 + 4)*x2)^2, ((x2^3/3 + (x1/2 + 2)*x2^2)*((1971*x1)/250 + (1971*x2)/250 + 3942/125))/(x2^2/2 + (x1 + 4)*x2) - (207*x1*x2^3*(x1 + x2 + x1*x2 + x2^2))/125 - (621*x1*x2^2*(x1*x2 + (x1*x2^2)/2 + x2^2/2 + x2^3/3))/125 + ((2*x2*(x1/2 + 2) + x2^2)*((1971*x2*(x1 + 4))/250 + (1971*x2^2)/500))/(x2^2/2 + (x1 + 4)*x2) - ((x2^3/3 + (x1/2 + 2)*x2^2)*((1971*x2*(x1 + 4))/250 + (1971*x2^2)/500)*(x1 + x2 + 4))/(x2^2/2 + (x1 + 4)*x2)^2, - x1 - 4/3]
Best regards
Stephan

Weitere Antworten (0)

Kategorien

Mehr zu Particle & Nuclear Physics finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by