[cig-commits] r6712 - geodyn/3D/MAG/trunk/src
wei at geodynamics.org
wei at geodynamics.org
Fri Apr 27 12:24:54 PDT 2007
Author: wei
Date: 2007-04-27 12:24:53 -0700 (Fri, 27 Apr 2007)
New Revision: 6712
Removed:
geodyn/3D/MAG/trunk/src/getgauss.f
geodyn/3D/MAG/trunk/src/unscramble.f
Log:
clean up, take out unscramble and getgauss subroutines
Deleted: geodyn/3D/MAG/trunk/src/getgauss.f
===================================================================
--- geodyn/3D/MAG/trunk/src/getgauss.f 2007-04-27 19:24:44 UTC (rev 6711)
+++ geodyn/3D/MAG/trunk/src/getgauss.f 2007-04-27 19:24:53 UTC (rev 6712)
@@ -1,82 +0,0 @@
- subroutine getgauss(l,m,alm,blm,glm,hlm,id)
-
-c this subroutine converts between Gauss coefficients (glm, hlm)
-c and the spherical harmonic coefficients (alm,blm)
-c of the magnetic potential at harmonic degree l
-c and order m used by the numerical dynamo model MAG.
-c
-c Inputs:
-c l,m: integer spherical harmonic degree and order;
-c
-c id: transformation direction:
-c id = 1 converts (alm,blm) to (glm,hlm)
-c id =-1 converts (glm,hlm) to (alm,blm)
-c
-c alm,blm: dimensionless spherical harmonic coefficients
-c of the radial magnetic field potential on the
-c outer boundary of the dynamo model. alm is the
-c coefficient of the real (cosine) part and blm
-c is for the imaginary (sine) part.
-c
-c glm,hlm: Gauss coefficients in nanotessla (nT).
-c glm multiplies with cos(m*phi)
-c hlm multiplies with sin(m*phi)
-c
-c The calculation uses SI units and nominal Earth values for
-c surface and core radii. The
-c conversion from dimensionless (alm,blm) to (glm,hlm) in
-c nT uses Earth values for density, rotation, electrical
-c conductivity, and uses the depth of the outer core as the
-c length scale.
-
- if (l .le. 0 .or. m .lt. 0 .or. m .gt. l) then
- write(6,'(''bad l or m in getgauss'')')
- return
- endif
-c
-c constants
-c
- pi=4.*atan(1.)
- pi4inv=1./(4.*pi)
- anano=1.0e+09
-c
-c nominal Earth parameters
-c
- sigma=6.0e+05 ! Core electrical conductivity in S/m
- rho=11.e+03 ! Core mean density in kg m^{-3}
- omega=7.2924e-05 ! Rotation angular frequency
- escale=sqrt(rho*omega/sigma) ! "Elsasser number" scaling
- re=6.371e+06 ! Earth radius in m
- rc=3.485e+06 ! Core radius in m
- ri=1.215e+06 ! Inner core radius in m
- cdepth=rc-ri ! Depth of outer core
- rcre=rc/re ! core/surface radius ratio
- rcsqinv=(cdepth/rc)**2 ! dim-less inverse squared core radius
-c
-c harmonic factors
-c
- fl=float(l)
- flp2=fl + 2.
- fl2p1=(2.*fl) + 1.
- fm=float(m)
- fact=sqrt(fl2p1*pi4inv)
- fact1=((-1.)**fm)*fl*fact
- if(m.gt.0) fact1=fact1*sqrt(2.)
- fact2=rcsqinv*(rcre**flp2)
-c
-c define MAG harmonic-Gauss harmonic conversion factors
-c conalm,conblm
-c
- conalm=fact1
- conblm=-fact1
-c
- if (id .gt. 0) then ! form Gauss coeffs in nT
- glm=anano*escale*fact2*conalm*alm
- hlm=anano*escale*fact2*conblm*blm
- return
- else ! form dimensionless fully normalized potential coeffs
- alm=glm/(anano*escale*fact2*conalm)
- blm=hlm/(anano*escale*fact2*conblm)
- return
- endif
- end
Deleted: geodyn/3D/MAG/trunk/src/unscramble.f
===================================================================
--- geodyn/3D/MAG/trunk/src/unscramble.f 2007-04-27 19:24:44 UTC (rev 6711)
+++ geodyn/3D/MAG/trunk/src/unscramble.f 2007-04-27 19:24:53 UTC (rev 6712)
@@ -1,60 +0,0 @@
- program unscramble
-c formulas in unscramble.f were taken from the prep.f and other mag files
-c it is used to unscramble the cc-files, converting the cc-file
-c output to arrays in harmonic degree l and order m.
-c It defines arrays al(lm) and am(lm) for lm=1,nlma,
-c which contain the l and m values we need to convert the cc-file contents
-
- integer lm, mc, lc
- real mclm(lm),la(lm),ma(lm)
- paramter (lmax, minc, nlma)
-
-
-c open cc file for reading
-
-c open (21,
-c input lmax, minc, and nlma from cc file header HERE:
- lmax=153
- minc=4
- nlma=32
-
-c define intermediate indices
- mmax=(lmax/minc)*minc
- nmaf=mmax+1
- nlaf=lmax+1
-
-c define the unscramble array mclm(lm)
- lm=0
- do 35 mc=1,nmaf,minc
- do 31 lc=mc,nlaf
- lm=lm+1
- mclm(lm)=mc
- 31 continue
- 35 continue
-
-c define la, ma arrays
- do 36 lm=1,nlma
- ma(lm)=mclm(lm)-1
-c define al in three terms
- tl1=lm+ma(lm)-1
- tl2=ma(lm)*(ma(lm)-minc)/(2.*minc)
- tl3=-ma(lm)*(lmax+1)/minc
- la(lm)=tl1+tl2+tl3
-c PRINT lm, la(lm), ma(lm) HERE
- 36 continue
-
-c if this works, then read the contents of the cc-file block in pairs,
-c and assign the new indices
-
-cc do 37 lm=1,nlma
-cc READ (cc-file) c1,c2
-c assign new indices
-cc l=la(lm)
-cc m=ma(lm)
-cc anew(l,m)=c1
-cc bnew(l,m)=c2
-cc 37 continue
-
-
-
- end
More information about the cig-commits
mailing list