Numerical solutions of equation

1 次查看（过去 30 天）
NILESH PANDEY
NILESH PANDEY 2018-7-22
评论： madhan ravi 2018-7-23
I'm trying to solve the following equation with respect to b
I need to find the value of b code is
clear all
E=1.5e8;
syms b
a0=-4e7/400;
T0=600;
T=300;
temp=1/(2*a0*(T0-T));
X=5050700286489334129256038400000/(4349164374701853503175341475971*b) + exp(erfinv(erf(pi/b)/2 + (633825300114114700748351602688*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b) - ((erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2)/2)^2)*exp(-(2*b*erfinv(erf(pi/b)/2 + (633825300114114700748351602688*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b) - ((erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2)/2) - b^2*1i)^2/(4*b^2))*(50507002864893341292560384/(4349164374701853503175341475971*b) + (50507002864893341292560384*(erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b)))/(4349164374701853503175341475971*b*((erf(pi/b) - (101014005729786682585120768*E + 4040560229191467303404830720000000)/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2))) + exp(erfinv(erf(pi/b)/2 + (633825300114114700748351602688*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b) - ((erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2)/2)^2)*exp(-(2*b*erfinv(erf(pi/b)/2 + (633825300114114700748351602688*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b) - ((erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2)/2) + b^2*1i)^2/(4*b^2))*(50507002864893341292560384/(4349164374701853503175341475971*b) + (50507002864893341292560384*(erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b)))/(4349164374701853503175341475971*b*((erf(pi/b) - (101014005729786682585120768*E + 4040560229191467303404830720000000)/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2))) + (5050700286489334129256038400000*(erf(pi/b) - (1267650600228229401496703205376*((5879788993030477*E)/73786976294838206464 + 459358515080506015625/144115188075855872))/(4349164374701853503175341475971*b)))/(4349164374701853503175341475971*b*((erf(pi/b) - (101014005729786682585120768*E + 4040560229191467303404830720000000)/(4349164374701853503175341475971*b))^2 + 1/5)^(1/2)) + 7737125245533627/9671406556917033397649408
vpasolve(X==temp,b)
but matlab(2016 version) gives the following error:
ans =
Empty sym: 0-by-1
equation also contain imaginary part which is written 1i in the expression of X
please help me for the same code for the calculation of X given below please have a look maybe this can help, code for the calculation of X is:
clear all
syms E b;
N0=1e6;
E0=1e8;
Eeff=-0.4e8;
a=19.3567;
p=0.1;
ed=0.4e-9;
efe=4.5*1e6*8.854e-12;
x=2*efe*(E-Eeff)/(a*b*p*sqrt(pi));
xm=erf(pi/b);
delta=0.05;
H_f=x-0.5*(x-xm+sqrt((x-xm).^2+4*delta));
theta=b*erfinv(H_f);
N1=N0*p;
J=N1*erf(pi/b);
p1=0.5*N1*exp(-b*b*0.5);
P_diel=2*ed*(E-E0);
p2=erf((2*theta+i*b*b)/(2*b))+erf((2*theta-i*b*b)/(2*b))+N1*erf(theta/b);%+P_diel;
P=p1+p2+P_diel;
X= diff(P,E);
  5 个评论
显示 3更早的评论
NILESH PANDEY
NILESH PANDEY 2018-7-23
Thanks for your time,I'll try some other way maybe series expansion. If I found the solution I'll share it.

请先登录，再进行评论。

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

采纳的回答

Walter Roberson
Walter Roberson 2018-7-23
%clear all
syms E b;
N0=1e6;
E0=1e8;
Eeff=-0.4e8;
a=19.3567;
p=0.1;
ed=0.4e-9;
efe=4.5*1e6*8.854e-12;
x=2*efe*(E-Eeff)/(a*b*p*sqrt(pi));
xm=erf(pi/b);
delta=0.05;
H_f=x-0.5*(x-xm+sqrt((x-xm).^2+4*delta));
theta=b*erfinv(H_f);
N1=N0*p;
J=N1*erf(pi/b);
p1=0.5*N1*exp(-b*b*0.5);
P_diel=2*ed*(E-E0);
p2=erf((2*theta+1i*b*b)/(2*b))+erf((2*theta-1i*b*b)/(2*b))+N1*erf(theta/b);%+P_diel;
P=p1+p2+P_diel;
X= diff(P,E);
Eval=1.5e8;
a0=-4e7/400;
T0=600;
T=300;
temp=1/(2*a0*(T0-T));
eqn = subs(X, E, Eval) - temp;
B = linspace(-10,10);
Y = double( subs(eqn, b, B) );
plot(B, Y)
vpasolve(eqn, b, [6 7])
answer is approximately 6.688061088729162561767129937816
  2 个评论
NILESH PANDEY
NILESH PANDEY 2018-7-23
I tried to run your code but Matlab gives the following error:
Error using symengine
DOUBLE cannot convert the input expression into a double array.
Error in sym/double (line 613)
Xstr = mupadmex('symobj::double', S.s, 0);
Error in sol (line 33)
Y = double( subs(eqn, b, B ) );
there is an error in line 33:
Y = double( subs(eqn, b, B ) );
is it because of Matlab version I'm using Matlab-2016?
NILESH PANDEY
NILESH PANDEY 2018-7-23
Thanks Sir for answer now I see above error is for the plotting, so I commented above lines in Matlab now matlab is able to solve above equation. Thanks again for your help and time!

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Conversion 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by