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]

```

• Prev by Date: Re: How do I create a parametric expression?
• Next by Date: a simple integral
• Previous by thread: Re: How do I create a parametric expression?
• Next by thread: a simple integral