12 views (last 30 days)

Show older comments

The coding are as below.

I got the error at LINE 5 where the erfc is not valid since there is a square root of complex number. But the calculation supposed to be like that, any other ways to write it without getting an error? I really appreciate some help here.

Coding:

y=linspace(0,6,150);

t=1; U=3; Pr=5; w=pi/3; B1=0.8; Gr=4.5; a1=1-Pr; B2=B1/a1; B3=Gr/(a1*B2);

H=1; r=1i*w*t; s1=B1+1i*w; s2=B1-1i*w; z1=1i*B2*Pr; z2=1i*B2*t;

u1=(1/4)*H*exp(imag(r)+y*sqrt(imag(s1))).*erfc(y/(2*sqrt(t))+sqrt(imag(s1))*t)+(1/4)*H*exp(imag(r)-y*sqrt(imag(s1))).*erfc(y/(2*sqrt(t))-sqrt(imag(s1)*t));

u2=(1/4)*H*exp(imag(-r)+y*sqrt(imag(s2))).*erfc(y/(2*sqrt(t))+sqrt(imag(s2))*t)+(1/4)*H*exp(imag(-r)-y*sqrt(imag(s2))).*erfc(y/(2*sqrt(t))-sqrt(imag(s2)*t));

u3=(B3/2)*exp(y*sqrt(B1))*erfc(y/(2*sqrt(t))+sqrt(B1*t))+(B3/2)*exp(-y*sqrt(B1))*erfc(y/(2*sqrt(t))-sqrt(B1*t));

u4=(1/2)*exp(-B2*t+y*sqrt(B1-B2))*erfc(y/(2*sqrt(t))-sqrt((B1-B2)*t))+(1/2)*exp(-B2*t-y*sqrt(B1-B2))*erfc(y/(2*sqrt(t))+sqrt((B1-B2)*t));

u5=B3*erfc(y/2*sqrt(Pr/t));

u6=(B3/2)*exp(-B2+y*imag(z1))*erfc((y/2)*sqrt(Pr/t)+imag(z2))+(B3/2)*exp(-B2-y*imag(z1))*erfc((y/2)*sqrt(Pr/t)-imag(z2));

u=u1+u2-u3+u4+u5-u6;

figure

plot(y,u,'blue')

title('Temperature profiles when Pr=0.7 and t=0.8')

xlabel('y')

ylabel('u')

David Goodmanson
on 22 Apr 2021

Hello Siti,

given the symmetry of u1 and u2 in the handwritten equations, it can't be the case that u1 has real argument for erfc and u2 has complex argument for erfc. It appears that they are all complex.

David Goodmanson
on 23 Apr 2021

Edited: David Goodmanson
on 5 May 2021

Hello Siti,

Although erfc does not take complex arguments as you know, the kummer U function is available and does take complex argument.

Kummer M and kummer U are partner functions where M is bounded at the origin, and U is not. In Matlab the kummer M function is denoted by hypergoem(a,b,z), and the kummer U function is denoted by kummerU. Matlab syntax is not very consistent here. Anyway, you can use

function y = erfcU(z)

% erfc from kummerU

%

% function y = erfcU(z)

yy = (1/sqrt(pi))*exp(-z.^2).*kummerU(1/2,1/2,z.^2);

a = real(z)>0 | (real(z)==0 & imag(z)>=0); % Re(z)>0, also positive y axis

y = (a.*yy + ~a.*(2-yy));

which works for complex argument.

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

Start Hunting!