Services & ResourcesWolfram ForumsMathGroup Archive

[Date Index] [Thread Index] [Author Index]

How to package an array generating code

• To: mathgroup at smc.vnet.net
• Subject: [mg68567] How to package an array generating code
• From: "Tonybony" <aneves at gmail.com>
• Date: Wed, 9 Aug 2006 23:57:07 -0400 (EDT)
• Sender: owner-wri-mathgroup at wolfram.com

Dear mathematica users,

I am trying to build my first package, read instructions on wolfram
website, but still the question remains. How to build a package from a
code that generates a "table" based on initial user input?

I have optimized the code for efficiecy, and do not see any faster
alternative for the code (attached below this email). The ideai of a
pakage was to share my work with other mathematica users and
portability.

I am happy to hear any suggestions and criticism on the code and a
solution for a package.

Thanks,
Antonio Neves

___________Mathematica Code____________
λ = 0.8; (*Wavelength in microns*)
xo = 1.5; (*Beam position in the x - axis in microns*)
yo = -1.5; (*Beam position in the x - axis in microns*)
zo = 1.5;(*Beam position in the x - axis in microns*)
NA = 1.25;(*Objective lens numerical aperture*)
f = 1700;(*Objective lens focal length in microns*)
Ï?b = 2500;(*Back objective lens apperture radius bin microns*)
nMedium = 1.4883;
nAir = 1.01;
k = (2*Pi*nMedium)/λ;
kÏ?o = k*Sqrt[xo^2 + yo^2];
kzo = k*zo;
kf = k*f;
fÏ? = f/Ï?b;
Ï?o = ArcTan[xo, yo];
Cosfio = Cos[Ï?o];
Sinfio = Sin[Ï?o];
Cosmax = Sqrt[1 - (NA/nMedium)^2];
Nmax = 30;

I1[(n_Integer)?Positive, m_Integer] :=
   NIntegrate[(Sqrt[t]*(m*LegendreP[n, m, t]*
         ((t*BesselJ[m + 1, kÏ?o*Sqrt[1 - t^2]])/
           Sqrt[1 - t^2] + (m*BesselJ[m,
             kÏ?o*Sqrt[1 - t^2]])/(kÏ?o*(t + 1))) -
        (m + n)*(-m + n + 1)*LegendreP[n, m - 1, t]*
         ((m*BesselJ[m, kÏ?o*Sqrt[1 - t^2]])/
           (kÏ?o*Sqrt[1 - t^2]) - BesselJ[m + 1,
           kÏ?o*Sqrt[1 - t^2]])))/(E^(fÏ?^2*t^2)*
       E^(I*k*t*zo)), {t, Cosmax, 1},
     AccuracyGoal -> 6, Compiled -> False,
     Method -> DoubleExponential, MaxRecursion -> 100]/
    E^fÏ?^2;
I2[(n_Integer)?Positive, m_Integer] :=
   NIntegrate[(Sqrt[t]*(LegendreP[n, m, t]*
         ((m*BesselJ[m, kÏ?o*Sqrt[1 - t^2]])/
           (kÏ?o*(1 + t)) - BesselJ[m + 1,
            kÏ?o*Sqrt[1 - t^2]]/Sqrt[1 - t^2]) -
        (n + m)*(n - m + 1)*LegendreP[n, m - 1, t]*
         (BesselJ[m, kÏ?o*Sqrt[1 - t^2]]/
          (kÏ?o*Sqrt[1 - t^2]))))/(E^(fÏ?^2*t^2)*
       E^(I*k*t*zo)), {t, Cosmax, 1},
     AccuracyGoal -> 6, Compiled -> False,
     Method -> DoubleExponential, MaxRecursion -> 100]/
    E^fÏ?^2;

GridAMIE = Table[{n, m}, {n, Nmax + 1}, {m, 0, n}];
MI1 = Apply[I1[#1, #2] & , GridAMIE, {2}];
MI2 = Apply[I2[#1, #2] & , GridAMIE, {2}];

GTM[(n_Integer)?Positive, (m_Integer)?NonNegative] :=
   4*Pi*I^(n - m)*Exp[(-I)*m*Ï?o]*I*kf*Exp[(-I)*kf]*
    Sqrt[nAir/nMedium]*
    Sqrt[((2*n + 1)/(4*Pi*n*(n + 1)))*
      ((n - m)!/(n + m)!)]*(MI1[[n,1 + Abs[m]]]*
      Cosfio + I*m*MI2[[n,1 + Abs[m]]]*Sinfio);
GTM[(n_Integer)?Positive, (m_Integer)?Negative] :=
   (-1)^m*Exp[-2*I*m*Ï?o]*GTE[n, -m];
GTE[(n_Integer)?Positive, (m_Integer)?NonNegative] :=
   -4*Pi*I^(n - m)*Exp[(-I)*m*Ï?o]*I*kf*Exp[(-I)*kf]*
    Sqrt[nAir/nMedium]*
    Sqrt[((2*n + 1)/(4*Pi*n*(n + 1)))*
      ((n - m)!/(n + m)!)]*(MI1[[n,1 + Abs[m]]]*
      Sinfio - I*m*MI2[[n,1 + Abs[m]]]*Cosfio);
GTE[(n_Integer)?Positive, (m_Integer)?Negative] :=
   (-1)^m*Exp[-2*I*m*Ï?o]*GTM[n, -m];
EAMIE[x_, y_, z_] := Module[{r, θ, Ï?, PartTE,
     PartTM, GridAMIE2, SumTE, SumTM},
    r = Sqrt[x^2 + y^2 + z^2];
     θ = ArcCos[z/Sqrt[x^2 + y^2 + z^2]];
     Ï? = ArcTan[x, y]; ParteTE[n_, m_] :=
      (-I)*GTE[n, m]*Sqrt[((2*n + 1)/(8*k*r*n*
           (n + 1)))*((n - m)!/(n + m)!)]*E^(I*m*Ï?)*
       BesselJ[n + 0.5, k*r]*
       {0, (-I)*m*(LegendreP[n, m, Cos[θ]]/
          Sin[θ]), (-(1 + n))*Cot[θ]*
          LegendreP[n, m, Cos[θ]] + (1 - m + n)*
          Csc[θ]*LegendreP[1 + n, m, Cos[θ]]};
     ParteTM[n_, m_] := (I/k)*GTM[n, m]*
       Sqrt[((2*n + 1)/(8*k*r*n*(n + 1)))*
         ((n - m)!/(n + m)!)]*E^(I*m*Ï?)*
       {I*n*(n + 1)*BesselJ[n + 0.5, k*r]*
         (LegendreP[n, m, Cos[θ]]/(k*r)),
        (-I)*0.5*(BesselJ[-0.5 + n, k*r] -
          BesselJ[0.5 + n, k*r]/(k*r) -
          BesselJ[1.5 + n, k*r])*((1 + n)*Cot[θ]*
           LegendreP[n, m, Cos[θ]] - (1 - m + n)*
           Csc[θ]*LegendreP[1 + n, m, Cos[θ]]),
        (-m)*0.5*(BesselJ[-0.5 + n, k*r] -
          BesselJ[0.5 + n, k*r]/(k*r) -
          BesselJ[1.5 + n, k*r])*(LegendreP[n, m,
           Cos[θ]]/Sin[θ])}; GridAMIE2 =
      Table[{n, m}, {n, Nmax}, {m, -n, n}];
     SumTE = Tr[Flatten[Apply[ParteTE[#1, #2] & ,
         GridAMIE2, {2}], 1], Plus, 1];
     SumTM = Tr[Flatten[Apply[ParteTM[#1, #2] & ,
         GridAMIE2, {2}], 1], Plus, 1];
     (SumTE + SumTM)*({x, y, z}/r)]; 

GTE[9, -5]
EAMIE[1.5, 1.5, 1.5]