using structures and iterators

*To*: mathgroup at yoda.physics.unc.edu*Subject*: using structures and iterators*From*: rclark at lpl.arizona.edu (Richard Clark x4971)*Date*: Thu, 7 Apr 94 07:47:33 MST

I am trying to find the propper use of Mma iterators and object structures to simplify a proceedure. I have two equivalent sets of equations and I need to find the relation between the sets of variables in which they are expressed to a desired order of aproximation. The variables being solved for are related to the extent that var1 is first order, var2 is 2nd order ... so that in the lower order partial solution the higher order varN (of either set) do not appear. x1,y1 ~m largest term in eqn1 is ~m x2,y2 ~m^2 largest term in eqn2 is ~m^2 x3,y3 ~m^3 largest term in eqn3 is ~m^3 and m << 1 I find the solutions to successively higher order: solve eqn1 for x1 to order m substitute this in eqn2 for and solve x2 to order m^2 substitute this in eqn1 for and solve x1 again, now to order m^2 substitute this in eqn3 for and solve x3 to order m^3 substitute this in eqn2 for and solve x2 again, now to order m^3 substitute this in eqn1 for and solve x1 again, now to order m^3 At present I am working only to 3rd order (which has already had a partial solution published). The following code works by performing each of the above steps explicitly. I want to extend this solution to higher order accuracy of solution for {x1,x2,x3,x4...}. Also this same type of successive aproximation solution will be used in several different steps in obtaining the complete solution to the response functions I am working on. Any suggestions would be apreciated. Richard Clark rclark at lpl.arizona.edu ------------------------ Now for the listings. (These variable names follow some traditional notation for the problem.) (e, Jj[1]) 1st order (Jj[1] is sometimes replaced with another 1st order (k, Jj[2]) 2nd order quantity 'ep' for reasons that were convienent (h, Jj[3]) 3rd order half a century ago) This first one has the advantage that it *does* work. I want to generalize it by using Mma's iterators instead of an editor! -------------------------------------------------------- (* hkslv: solve (ep,Jn) <--> (e,k,h) *) (* get the 2 equivalent sets of equations *) <<rslt/ra.res <<rslt/rahk.res (* define the order of the terms *) qhrul={q->aa q,ep->aa ep,e->aa e,k->aa^2 k,h->aa^3 h,Jj[n_]->aa^n Jj[n]} ehkrul={ee[n_]->aa^n ee[n], kk[n_]->aa^(n+1) kk[n], hh[n_]->aa^(n+2) hh[n]} (* after substitution express e,k,h as sums of various orders *) ehksub={e->ees, k->kks, h->hhs} ees=Sum[ee[n],{n,1,maxord}] kks=Sum[kk[n],{n,1,maxord-1}] hhs=Sum[hh[n],{n,1,maxord-2}] (* make copies of the 2 systems containing terms up to a desired order *) Subser[ord_] := Block[{i}, Do[sshktmp[i]=Normal[Series[sshk[i]/.ehksub/.qhrul/.ehkrul,\ {aa,0,ord}]/.aa->1],{i,1,maxord}]; Do[sstmp[i]=Normal[Series[ss[i]/.qhrul,\ {aa,0,ord}]/.aa->1],{i,1,maxord}]; sstm=Table[sstmp[i]==sshktmp[i],{i,1,maxord}]; ] (* 1st order *) Subser[1] e1=Solve[sstm[[1]],ee[1]] ee[1]=ee[1]/.e1[[1]] (* 2nd order *) Subser[2] k1=Solve[sstm[[2]],kk[1]] kk[1]=kk[1]/.k1[[1]] e2=Solve[sstm[[1]],ee[2]] ee[2]=ee[2]/.e2[[1]] (* 3rd order *) Subser[3] h1=Solve[sstm[[3]],hh[1]] hh[1]=hh[1]/.h1[[1]] k2=Solve[sstm[[2]],kk[2]] kk[2]=kk[2]/.k2[[1]] e3=Solve[sstm[[1]],ee[3]] ee[3]=ee[3]/.e3[[1]] ees=Expand[ees] kks=Expand[kks] hhs=Expand[hhs] Save["hkslv.res", ees, kks, hhs] -------------------------------------------------------- This is what I've come up with so far. It looks like it should be close to what I want but doesn't work. (* hkslv: solve (ep,Jn) <--> (e,k,h) (21-30) *) <<rslt/ra.res <<rslt/rahk.res maxord=3 (* define the order of terms *) qhrul={q->aa q,ep->aa ep,e->aa e,k->aa^2 k,h->aa^3 h,Jj[n_]->aa^n Jj[n]} (* make a list of the variables to solve for *) ehk={e, k, h} (* called as sln=Solser[3, sshk, ss, ehk] *) (* set1 is the system expressed in varset *) Solser[ord_,set1_,set2_,varset_]:=Block[{}, vsub=Table[varset[[i]]->v[i],{i,1,maxord}]; vrul={vn[n_,m_]->aa^m vn[n,m]}; (* ideally the solution will end up in v *) (* I count the order of the vn differently from [ehk]s in version 1 so that the second subscript is actually the order of the term *) Do[v[i]=Sum[vn[i,j],{j,i,maxord}],{i,1,maxord}] Do[Block[{}, (* these 3 statements are the equivalent of Subser in version 1 *) Do[s1s[i]=Normal[Series[set1[i]/.vsub/.qhrul/.vrul,\ {aa,0,ord}]/.aa->1],{i,1,iord}]; Do[s2s[i]=Normal[Series[set2[i]/.qhrul,\ {aa,0,ord}]/.aa->1],{i,1,iord}]; stm=Table[s1s[i]==s2s[i],{i,1,iord}]; Do[Block[{}, slv=Solve[stm[[iord-jord+1]],vn[iord-jord+1,iord]]; vn[iord-jord+1,iord]=vn[iord-jord+1,iord]/.slv[[1]]; ],{jord,1,iord}]; ],{iord,1,ord}] ] -------------------------------------------------------- Here is ./rslt/ra.res maxord = 3 ss/: ss[0] = (-4*ep^2)/45 - (512*ep^3)/2835 + (4*ep^2*q)/35 ss/: ss[1] = (-2*ep)/3 - (44*ep^2)/63 - (1286*ep^3)/945 + (3*ep*q)/7 + (25*ep^2*q)/21 - (4*ep*q^2)/21 - (9*ep*Jj[2])/14 - (25*q*Jj[2])/168 ss/: ss[2] = (-16*ep^2)/35 - (1376*ep^3)/1155 + (4*ep*q)/7 + (236*ep^2*q)/231 - (18*ep*q^2)/77 - Jj[2] - (548*ep*Jj[2])/231 - (18*q*Jj[2])/77 ss/: ss[3] = (-16*ep^3)/33 + (160*ep^2*q)/231 - (100*ep*q^2)/231 - (20*ep*Jj[2])/11 + (25*q*Jj[2])/33 - Jj[3] -------------------------------------------------------- and ./rslt/rahk.res maxord = 3 sshk/: sshk[0] = (-4*e^2)/45 - (52*e^3)/567 - (32*e*k)/315 sshk/: sshk[1] = (-2*e)/3 - (23*e^2)/63 - (4*e^3)/27 + (2*h)/21 - (8*k)/21 - (152*e*k)/315 sshk/: sshk[2] = (12*e^2)/35 + (4*e^3)/11 + (192*h)/385 + (32*k)/35 + (32*e*k)/105 sshk/: sshk[3] = (-40*e^3)/231 - (80*h)/231