diff options
Diffstat (limited to 'maxima_code/w2wl.max')
| -rw-r--r-- | maxima_code/w2wl.max | 288 |
1 files changed, 288 insertions, 0 deletions
diff --git a/maxima_code/w2wl.max b/maxima_code/w2wl.max new file mode 100644 index 0000000..158c900 --- /dev/null +++ b/maxima_code/w2wl.max @@ -0,0 +1,288 @@ +writefile("session.out")$ + +dist(z):=(MATRIX([1,z],[0,1]))$ +lens(f):=(MATRIX([1,0],[-1/f,1]))$ + +abcd(optics_system):=block([tmp:MATRIX([1,0],[0,1])], (for i:1 thru length(optics_system) do tmp:optics_system[i].tmp ), return(tmp))$ + +prop(q0,abcd):=(abcd[1][1]*q0+abcd[1][2])/(abcd[2][1]*q0+abcd[2][2])$ + +wr2q(w,r):=1/(1/r-%I*lamb/%PI/(w^2))$ +q2w(q):=sqrt(-lamb/%PI/imagpart(1/q))$ +q2r(q):=block([rp:(realpart(1/q))], + (if (rp = 0) + then + (print("#dbg: R=1/0, will substitute R with something big"), + return(1e100)) + else + return(1/rp) + ) +) $ + + +/* find q0 for such q */ +q2w0(q):=block([w,r,w0], + w:q2w(q),r:q2r(q), + w0:(w/sqrt(1+(%Pi*w^2/lamb/r)^2)), + return (w0) + )$ +q2z0(q):=block([w,r], + w:q2w(q),r:q2r(q), return(r/(1+(lamb*r/%Pi/w^2)^2)) + )$ +q2q0(q):=%I*%Pi*q2w0(q)^2/lamb$ + +q0f2z_for_colim_lens(q0,f):=f+%Pi*q2w(q0)^2/lamb$ + +qf2colim_lens_posit(qstart,f):= block([q0,z0,zd], + q0:q2q0(qstart), + z0:q2z0(qstart), + zd:q0f2z_for_colim_lens(q0,f), + return(zd-z0) +)$ + + +check_real(x):=if featurep(x[1],real) then if featurep(x[2],real) then if featurep(x[3],real) then true else false $ + +check_posit(x):=if (part(x[1],2)>0) then if (part(x[2],2)>0) then if (part(x[3],2)>0) then true else false $ + +check_sm1(x):=if (part(x[1],2)<ltotal) then if (part(x[2],2)<ltotal) then if (part(x[3],2)<ltotal) then true else false $ + +select_realsolutions(x):=block([temp:[]], for i:1 thru length(x) do (if check_real(x[i]) then temp:float(append(temp,[x[i]])) ), return(temp) ) $ + +select_posit_solutions(x):=block([temp:[]], for i:1 thru length(x) do (if check_posit(x[i]) then temp:float(append(temp,[x[i]])) ), return(temp) ) $ + +select_sm1_solutions(x):=block([temp:[]], for i:1 thru length(x) do (if check_sm1(x[i]) then temp:float(append(temp,[x[i]])) ), return(temp) ) $ + +/*=========================================================*/ +/* + match one Gauss beam parameter q0(w0,r0) + to another one qf(wf,qf) + with the help of two known lens. + distance between q0 and qf is 'ltotal' + q0>|----z1-----(f1)------z2-------(f2)------z3--------|>qf + ^ ^ + |-------------------ltotal-------------------------| + the system should looks like above. + */ +/*=========================================================*/ +qffql2z1z2z3(q0,f1,f2,qf,ltotal):=block( + [ + m1,m2,m3,m4,m5,z1,z2,z3, + optics, mfull, + q0r,qoi, qfr,qfi, + eqt,eq1,solut, + eq1,eq2,eq3, + eqn1, eqn2,eqn3, + eqnev1, eqnev2, eqnev3, + slt, rsol + ], + + m1:dist(z1), + m2:lens(f1), + m3:dist(z2), + m4:lens(f2), + m5:dist(z3), + + optics:[m1,m2, m3, m4, m5], + mfull:abcd(optics), + + eqt:qf=prop(q0,mfull), + + q0r:float(realpart(q0)), + q0i:float(imagpart(q0)), + qfr:float(realpart(qf)), + qfi:float(imagpart(qf)), + + eq1:ev(eqt, q0=q0r+q0i*%I, qf=qfr+qfi*%I), + + solut:solve([eq1], [z3]), + + eqn1:imagpart(solut[1]), + eqn2:realpart(solut[1]), + eqn3:z1+z2+z3=ltotal, + + eqnev1:ev(eqn1,numer), + eqnev2:ev(eqn2,numer), + eqnev3:ev(eqn3,numer), + + + slt:solve([eqnev1,eqnev2,eqnev3],[z1,z2,z3]), + + print("slt=",slt), + + rsol:select_realsolutions(slt), + rsol:select_posit_solutions(rsol), + rsol:select_sm1_solutions(rsol), + + print("rsol=",rsol), + + + return(rsol) +); +/*=========================================================*/ +/* + Zero ewaste to zero waste transformatin with lense + q0|------------(f)-----------|qf0 + ^ ^ ^ + |----z--------| | + | | + |-------------l------------| + we are lookking for z1,l +*/ + q0fq02z1(q0,f,qf):=block([f0,z1,w0tmp:q2w(q0),wftmp:q2w(qf)], + f0:%Pi*w0tmp*wftmp/lamb, + z1:f+w0tmp/wftmp*sqrt(f^2-f0^2), + return(z1) + )$ + q0fq02z2(q0,f,qf):=block([f0,z2,w0tmp:q2w(q0),wftmp:q2w(qf)], + f0:%Pi*w0tmp*wftmp/lamb, + z2:f+wftmp/w0tmp*sqrt(f^2-f0^2), + return(z2) + )$ + q0fq02l(q0,f,qf):=q0fq02z1(q0,f,qf)+q0fq02z2(q0,f,qf)$ +/*=========================================================*/ +/* + for back propogation we need to change R -> -R +*/ +q2bq(q):=wr2q(q2w(q),-q2r(q)); +/*=========================================================*/ +/*=========================================================*/ +/*=========================================================*/ +/*=========================================================*/ + + + +/*=========================================================*/ +/* lets solve final waste to parrallel beam, + (actually back propogation) +*/ +/*=========================================================*/ + +optics_after_parral_beam:[lens(-.111),dist(.05)]$ + +ltotal:4.85$ +lamb:1.064e-6$ +/* modecleaner waste */ +w0:25.63e-6$ +r0:1e100$ + q0:ev(wr2q(w0,r0),numer)$ + +q:ev(prop(q0,abcd([dist(.05),lens(-.111)])),numer); +dist2waste:ev(-q2z0(q),numer); + +quit(); + + +/* OPA waste */ +wend:37.09e-6$ +rend:1e100$ + qend:wr2q(wend,rend)$ +/* +qf:prop(qend,abcd(reverse(optics_after_parral_beam)))$ +float(realpart(qf))$ +float(imagpart(qf))$ +*/ +qf:qend$ +f_col:.203$ + +/*finally position of collimated lens */ +z_col:float(qf2colim_lens_posit(qf,f_col))$ +if (z_col <0) + then (print("it is not possible to get collimated beam z_col=",z_col), + quit() +) $ +qf:prop(qf,abcd([dist(z_col),lens(f_col)]))$ + + + +f1:.075$ +f2:.075$ + +dist_betw_f2_f_col:.40$ + +qf:prop(qf,abcd([dist(dist_betw_f2_f_col),lens(f2)]) )$ + +/* below something wrong */ +dist2waste:ev(-q2z0(qf),numer)$ +rad_now:ev(q2r(qf),numer)$ +qf:ev(q2q0(qf),numer)$ +/*===========================================*/ + + +space_between_wastes:ev(q0fq02l(q0,f1,qf), numer)$ +dist2first_lens:ev(q0fq02z1(q0,f1,qf),numer)$ +dist2first_lens_tmp:ev(q0fq02z2(q0,f1,qf),numer)$ + +z1:ev(dist2first_lens,numer)$ +z2:ev(space_between_wastes+dist2waste,numer)$ +z3:ev(z2+dist_betw_f2_f_col,numer)$ + +ltotal:z3+z_col$ + +/* +qf:wr2q( -q2r(qf), q2w(qf) )$ +*/ + + +/* +rsol:qffql2z1z2z3(q0,f1,f2,qf,ltotal-z_col)$ + print("rsol=",rsol), +*/ + + + with_stdout("answ.txt", + print("##########################################"), + print( "lambda=",lamb), + print( "Ltot=", ltotal), + print( "r0=", r0), + print( "w0=", w0), + print(), + print( "lns1.abcd=abcd_lens(", f1 , ")" ), + print( "lns1.x=", z1 ), + print( "lns2.abcd=abcd_lens(", f2 , ")" ), + print( "lns2.x=", (z2) ), + print( "lns3.abcd=abcd_lens(", f_col , ")" ), + print( "lns3.x=", (z3) ), + print("wf=", wend), + print("rf=", rend), + /* + qfinal:prop(q0,ev(mfull)), + wfinal:float(q2w(qfinal)), + rfinal:float(q2r(qfinal)), + print("#final waste=", wfinal), + print("#I was aiming to final waste"), + print("#final radius=", rfinal), + print("#I was aiming to final radius"), + */ + print("##########################################") + )$ + system("cat answ.txt")$ + system("xterm -e octave abcd.m") $ + +quit()$ + +/* +solve([eqn1,eqn2,eqn3],[z1,z2,z3]); +*/ + + + +/* + + +z3:ltotal-z1-z2$ +fabcd:ev(mfull)$ +f(z1,z2):=ev(float(abs(qf)-abs(prop(q0,fabcd))))$ +fr(z1,z2):=abs(f(z1,z2)); +fi(z1,z2):=imagpart(f(z1,z2)); +fr(.2,.4); +fi(.2,.4); + +plot3d(f(z1,z2), [z1,0,1],[z2,0,1], ['colour_z,true], [plot_format,gnuplot], +[run_viewer,false]); + +printreal(x):=if featurep(x[1],real) then if featurep(x[2],real) then if featurep(x[3],real) then x; + +printreal(slt[4]);float(ev(realpart(solut),number)); +float(ev(imagpart(solut),number)); +*/ |