adel-00

95 Reputation

9 Badges

9 years, 303 days

MaplePrimes Activity


These are replies submitted by adel-00

@acer 

Thanks for your responds

first: y is a complex; lets say y=a+ib, So y^2 is Re(y)^2+Im(y)^2

the isuue is all integrations contains deltao(unknown), so how can I plot y^2 against deltao 

is there  a way to do it without loop?

 

it is not at each point of the vectror yr is assigned to sls

for example for 0.3 in vertical axis should be verticle coulmn (200 )

@Carl Love 

Thanks

But why in matlab produce the figure different?

what shall i add in maple in the  plot option to get the same result..

 

@Ramakrishnan 

Thanks

xa and xb are functions of t

P1:=logplot([xa,xb],0..1,axes=boxed,title=tit,color=black,font=[2,3,18],thickness=2,tickmarks=[3,2],titlefont=[SYMBOL,14],font=[1,1,18],linestyle=[1,3,2]);

but still the same issue not like the matlab plot

@acer 

Hi 

I found error in the integration so i did correct the expressions but I got computation interrupted:

 

here the code below plz I appreciated if you can give me comments

restart: # d is real (plot of Spec agianst d)
    
kernelopts(version):    


term1:=(exp((1+I*d)*Gamma*tau)*(1-cos(2*Omega1))+cos(2*Omega1)*exp((1+I*d)*Gamma*t0)-exp(-(1+I*d)*Gamma*t0))/(2*(1+I*d)):
term2:=-1/2*evalf(Int(exp((1+I*d)*Gamma*x)*cos(Omega1*(1+erf(x))),x=-3..3)):
J1:=evalf(term1+term2):

J1mod:=(Re(J1))^2+(Im(J1))^2:

###### J2#########################
A2:=Int(exp((1+I*d)*Gamma*x)*sin(Omega1*(1+erf(x))),x=-3..3):
A3:=sin(2*Omega1)*(exp((1+I*d)*Gamma*tau)-exp((1+I*d)*Gamma*t0))/(1+I*d):
J2:=-I*(A2+A3):
J2mod:=(Re(J2))^2+(Im(J2))^2:
#J3 same as J1differ in sign
term1:=(exp((1+I*d)*Gamma*tau)*(1+cos(2*Omega1))-cos(2*Omega1)*exp((1+I*d)*Gamma*t0)+exp(-(1+I*d)*Gamma*t0))/(2*(1+I*d)):
term2:=0.5*int(exp((1+I*d)*Gamma*x)*cos(Omega1*(1+erf(x))),x=-3..3):
J3:=term1+term2:
J3mod:=(Re(J3))^2+(Im(J3))^2:
J4:=-J2:
J4mod:=(Re(J4))^2+(Im(J4))^2:

 


#Spec:=J1mod+J2mod:
Spec:=J1mod*cos(theta/2)^2+J2mod+J3mod*sin(theta/2)^2 -0.5*Re(J3*J4*sin(theta)*exp(I*phi))+0.5*Re(J1*J4*sin(theta)*exp(-I*phi)):

with(plots):
tau:=Pi:tau0:=1:Omega1:=10;Gamma:=1:t0:=3:theta:=Pi/3:phi:=Pi/3:
                          Omega1 := 10
tit:=sprintf("W=%g,G=%g,t=%g,q=%g,f=%g",Omega1,Gamma,tau,theta,phi);
         tit := "W=10,G=1,t=3.14159,q=1.0472,f=1.0472"
P1:=plot(Spec,d=-350..350,axes=boxed,numpoints=200,title=tit,color=black,
         font=[2,3,18],thickness=2,tickmarks=[3,3],gridlines=false,
         labels=["",""],
         titlefont=[SYMBOL,14],font=[1,1,18],linestyle=1);
Warning,  computation interrupted
Normalize:= proc(P::specfunc(anything, PLOT))
  local A,Smax1;
  A:= op([1,1], P);
  Smax1:= max(A[..,2]);
  if A::list then A:= Matrix(A) end if;
  A[..,2]:= A[..,2]/Smax1;
  subsop([1,1]= A, P);
end proc:
P1:= Normalize(P1):
for kk from 2 to 3 do
  tau:= 0.2*kk*Pi;
  P||kk:= plot(Spec,d=-350..350,axes=boxed,numpoints=200,title=tit,color=black,
               font=[2,3,18],thickness=2,tickmarks=[3,3],gridlines=false,
               labels=["",""],
               titlefont=[SYMBOL,14],font=[1,1,18],linestyle=1);
  P||kk:= plottools:-translate(Normalize(P||kk), 0, kk-1)
od:
display([P||(1..3)],view=[-10..10,0..3]);

@acer 

it is really run quicker many thanks

@acer 

              "Normalised spectrum for all cases"
restart:
------------------------- Defining the nature of the variables used ----------------------
assume(theta,real):assume(phi,real):assume(d,real):assume(tau,real):


                           tau0 := 1
J1
term1:=(exp((1+I*d)*Gamma*tau)*(1-cos(2*Omega1))+cos(2*Omega)*exp((1+I*d)*Gamma*t0)-exp(-(1+I*d)*Gamma*t0))/(2*(1+I*d)):

term2:=-0.5*evalf(int(exp((1+I*d)*Gamma*x)*cos(Omega1*(1+erf(tau))),x=-1..1)):
J1:=(term1+term2):
J1mod:=(Re(J1))^2+(Im(J1))^2:

###### J2#########################
A2:=evalf(int(exp((1+I*d)*Gamma*x)*sin(Omega1*(1+erf(tau))),x=-1..1)):
A3:=sin(2*Omega1)*(exp((1+I*d)*Gamma*tau)-exp((1+I*d)*Gamma*t0))/(1+I*d):

J2:=-I*(A2+A3):
######################

J2mod:=(Re(J2))^2+(Im(J2))^2:

J3 same as J1differ in sign
term1:=(exp((1+I*d)*Gamma*tau)*(1+cos(2*Omega1))-cos(2*Omega)*exp((1+I*d)*Gamma*t0)+exp(-(1+I*d)*Gamma*t0))/(2*(1+I*d)):

term2:=0.5*evalf(int(exp((1+I*d)*Gamma*x)*cos(Omega1*(1+erf(tau))),x=-1..1)):

 

J3:=term1+term2:
J3mod:=(Re(J3))^2+(Im(J3))^2:
J4 =- J2 
J4:=-J2:
J4mod:=(Re(J4))^2+(Im(J4))^2:
calculate the spectrum
Spec:= unapply((J1mod+J2mod,d)):
#Spec:= unapply((J1mod*cos(theta/2)^2+J2mod+J3mod*sin(theta/2)^2-0.5*Re(J3*J4*sin(theta)*exp(I*phi))+0.5*Re(J1*J4*sin(theta)*exp(-I*phi))),d):
with(plots):

Omega:=1:tau:=Pi:tau0:=1:Omega1:=0.5*Omega*tau0*sqrt(Pi):Gamma:=0.05:t0:=0.75:theta:=0:phi:=0:
#tit:=sprintf("l=%g,W=%g,G=%g,q=%g,f=%g",lambda,Omega,Gamma,theta,Phi):
P1:=plot((Spec),-350..350,axes=boxed,numpoints=200,title=tit,color=black,font=[2,3,18],thickness=2,tickmarks=[3,3],titlefont=[SYMBOL,14],font=[1,1,18],linestyle=1);

Normalize:= proc(P::specfunc(anything, PLOT))
local A,Smax1;
A:= op([1,1], P);
Smax1:= max(A[..,2]);
if A::list then A:= Matrix(A) end if;
A[..,2]:= A[..,2]/Smax1;
subsop([1,1]= A, P)
end proc:

P1:= Normalize(P1):
for kk from 2 to 5 do
tau0:= kk*Pi;
P||kk:= plot(Spec,-350..350,axes=boxed,numpoints=200,title=tit,color=black,font=[2,3,18],thickness=2,tickmarks=[3,3],titlefont=[SYMBOL,14],font=[1,1,18],linestyle=1);
P||kk:= plottools:-translate(Normalize(P||kk), 0, kk-1)
od:

display([P||(1..5)],view=[-350..350,0..5]);

that is the code that the integration in the program.. but it isn;t work

 

@acer 

hi if i make changes where

evalf(int(exp((5+I*d)*x)*sin(1+erf(x)),x=-1.2..1.2));

as d is areal number

in my opinoin I try to:

sin(erf) is equal to 0.5( exp(erf)-iexp(-erf))where erf(x)= 2/sqrt(pi) sum(x^2n+1)n!(2n+1)

so if I find the term exp(ax+bx^3) will be brilliant 

 

@Christian Wolinski 

to find integration of sin(erf(t)) e^(a+ib)t I need explicit formula

 

many thanks of ur response

 

@Carl Love I am not sure if I can exapnd it either with Bessel function or somthing simmilar

if i integrate sin(erf(t)) e^(a+ib)t

@acer thanks I got it

@Preben Alsholm 

and with that

restart:Digits:=70:
------------------------- Defining the nature of the variables used ----------------------
N:=0:M:=0:N1:=1+N:w:=10:


ini1:= x(0)=0.5,y(0)=0.5,z(0)=0;
            ini1 := x(0) = 0.5, y(0) = 0.5, z(0) = 0
var:={x(t),y(t),z(t)}: 
dsys:={diff(z(t),t)=-(N1+M*cos(2*w*t))*z(t)-1+f*(x(t)+y(t)), diff(x(t),t)=-(N1-I*w-2*M*exp(-2*I*w*t))*x(t)-f*(N1+(z(t)))-2*f*M*exp(2*I*w*t),diff(y(t),t)=-(N1+I*w-2*M*exp(2*I*w*t))*y(t)-f*(N1+(z(t)))-2*f*M*exp(-2*I*w*t)}:
zd:=subs(dsys,diff(z(t),t));
               zd := -z(t) - 1 + f (x(t) + y(t))
res:=dsolve(dsys union {x(0)=0.5,y(0)=0.5,z(0)=0},numeric,output=listprocedure):
Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)
#tit:=sprintf("F=%g,N=%g",f,N):
plot3d(z(t,f), f=-1..1, t=0..3, grid=[29,49]);

@acer 

Many thanks for ur response 

If I want to plot z(t) 

CodeTools:-Usage( plot3d(zfun(t,f), f=-1..1, t=0..3, grid=[29,49]) );

@Preben Alsholm 

N:=0:M:=0:N1:=1+N:w:=10:


ini1:= x(0)=0.5,y(0)=0.5,z(0)=0;
            ini1 := x(0) = 0.5, y(0) = 0.5, z(0) = 0
var:={x(t),y(t),z(t)}: 
dsys:={diff(z(t),t)=-(N1+M*cos(2*w*t))*z(t)-1+f*(x(t)+y(t)), diff(x(t),t)=-(N1-I*w-2*M*exp(-2*I*w*t))*x(t)-f*(N1+(z(t)))-2*f*M*exp(2*I*w*t),diff(y(t),t)=-(N1+I*w-2*M*exp(2*I*w*t))*y(t)-f*(N1+(z(t)))-2*f*M*exp(-2*I*w*t)}:
zd:=subs(dsys,diff(z(t),t));
                 zd := -z(t) - 1 + x(t) + y(t)
res:=dsolve(dsys union {x(0)=0.5,y(0)=0.5,z(0)=0},numeric,output=listprocedure):

 

I mean plot 3d of z(t),t,f 

1 2 3 4 5 6 7 Last Page 1 of 9