Scrambler  1
clumps.f90
Go to the documentation of this file.
00001 !#########################################################################
00002 !               
00003 !    Copyright (C) 2003-2012 Department of Physics and Astronomy,
00004 !                            University of Rochester,
00005 !                            Rochester, NY
00006 !
00007 !    clumps.f90 is part of AstroBEAR.
00008 !
00009 !    AstroBEAR is free software: you can redistribute it and/or modify    
00010 !    it under the terms of the GNU General Public License as published by 
00011 !    the Free Software Foundation, either version 3 of the License, or    
00012 !    (at your option) any later version.
00013 !
00014 !    AstroBEAR is distributed in the hope that it will be useful, 
00015 !    but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 !    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00017 !    GNU General Public License for more details.
00018 !
00019 !    You should have received a copy of the GNU General Public License
00020 !    along with AstroBEAR.  If not, see <http://www.gnu.org/licenses/>.
00021 !
00022 !#########################################################################
00025 
00029 
00032 
00036 
00039 MODULE Clumps
00040    USE GlobalDeclarations
00041    USE DataDeclarations
00042    USE PhysicsDeclarations
00043    USE CommonFunctions
00044    USE Shapes
00045    USE Perturbation
00046    USE VectorPerturbation
00047    USE DataInfoOps
00048    USE EOS
00049    USE BE_MODULE
00050    USE ObjectDeclarations
00051    IMPLICIT NONE
00053    TYPE ClumpDef
00054       REAL(KIND=qPREC) :: density=10   !Clump peak density
00055       REAL(KIND=qPREC) :: pressure=0   !Clump Pressure
00056       REAL(KIND=qPREC) :: temperature=.1      !Clump temperature
00057       REAL(KIND=qPREC), DIMENSION(3) :: velocity= (/0,0,0/)  !Clump velocity 
00058       REAL(KIND=qPREC), DIMENSION(3) :: position = (/0,0,0/)  !Clump velocity 
00059       REAL(KIND=qPREC) :: radius=1
00060       REAL(KIND=qPREC) :: phi=0
00061       REAL(KIND=qPREC) :: theta=0
00062       REAL(KIND=qPREC) :: thickness=.1 !thickness of clump smoothing region
00063       !      REAL(KIND=qPREC) :: B_tor=0     !Maximum Bfield for toroidal configuration
00064       !      REAL(KIND=qPREC) :: B_pol=0     !Maximum Bfield for poloidal configuration
00065       REAL(KIND=qPREC) :: B_phi=0
00066       REAL(KIND=qPREC) :: B_theta=0
00067       REAL(KIND=qPREC) :: B_toroidal=0     !Maximum Bfield for toroidal configuration
00068       REAL(KIND=qPREC) :: B_poloidal=0     !Maximum Bfield for poloidal configuration
00069       REAL(KIND=qPREC) :: omega=0     !Angular velocity for solid body rotation
00070       INTEGER :: RefineLevel=MaxDepth
00071       INTEGER :: SubSample=1
00072       INTEGER :: iTracer=0                                !Clump Tracer
00073       INTEGER :: density_profile=0
00074       REAL(KIND=qPREC) :: m2A=0 !Amplitude of an m=2 density perturbation [  rho = rho_0*(1+m2A*cos(2*phi))   ]
00075       LOGICAL :: PersistInBoundaries(3,2) = .false.
00076       INTEGER :: id
00077       TYPE(ShapeDef), POINTER :: Shape => Null()
00078       TYPE(PerturbationDef), POINTER :: DensityPerturbation => Null()
00079       TYPE(VectorPerturbationDef), POINTER :: VelocityPerturbation => Null()
00080       TYPE(PerturbationDef), DIMENSION(:), POINTER :: MagneticPerturbation => Null()
00081       INTEGER :: ObjId
00082    END TYPE ClumpDef
00083 
00084    !new declaration
00085    TYPE pClumpDef
00086       TYPE(ClumpDef), POINTER :: ptr
00087    END TYPE pClumpDef
00088    TYPE(pClumpDef) :: pClump   
00089    !
00090 
00091    INTEGER, PARAMETER :: UNIFORM = 0, BE_PROFILE = 1
00092    SAVE
00093 CONTAINS
00094 
00095 
00098    SUBROUTINE CreateClump(Clump, density, pressure, position)
00099       TYPE(ClumpDef), POINTER :: Clump
00100       REAL(KIND=qPREC), OPTIONAL :: density, pressure, position(3)
00101       ALLOCATE(Clump) 
00102       ALLOCATE(clump%Shape)
00103       IF (Present(density)) Clump%density=density
00104       IF (Present(pressure)) Clump%pressure=pressure
00105       IF (Present(position)) Clump%position=position
00106       CALL AddClumpToList(Clump)
00107       CALL UpdateClump(Clump)
00108    END SUBROUTINE CreateClump
00109 
00112    SUBROUTINE UpdateClump(Clump)
00113       TYPE(ClumpDef), POINTER :: Clump
00114       CALL SetShapeType(Clump%Shape, SPHERE, (/Clump%radius,0d0,0d0/))
00115       CALL SetShapeOrientation(Clump%Shape, 0d0, Clump%theta, Clump%phi)
00116       Clump%Shape%Position=Clump%position
00117       Clump%Shape%velocity=Clump%velocity
00118       Clump%Shape%t0=levels(0)%tnow
00119       CALL SetShapeBounds(Clump%Shape)         
00120    END SUBROUTINE UpdateClump
00121 
00122    SUBROUTINE AddClumpToList(Clump)
00123       TYPE(ClumpDef), POINTER :: Clump
00124       TYPE(ObjectDef), POINTER :: Object
00125       Clump%ObjId = ObjectListAdd(Object,CLUMPOBJ)
00126       pClump%ptr => Clump
00127       len = size(transfer(pClump, dummy_char))
00128       ALLOCATE(Object%storage(len))
00129       Object%storage = transfer(pClump, Object%storage)
00130    END SUBROUTINE AddClumpToList
00131 
00132    SUBROUTINE DestroyClump(Clump)
00133       TYPE(ClumpDef), POINTER :: Clump
00134       TYPE(ObjectDef),POINTER :: Object
00135       INTEGER :: id
00136       id=Clump%ObjId
00137       Object => ObjectListFind(id)
00138       IF (ASSOCIATED(Object) .AND. Object%type == CLUMPOBJ) THEN
00139          pClump = transfer(Object%storage,pClump)
00140          DEALLOCATE(pClump%ptr)
00141          NULLIFY(pClump%ptr)
00142          CALL ObjectListRemove(id)
00143       ENDIF
00144       NULLIFY(Clump)
00145    END SUBROUTINE DestroyClump
00146 
00147    SUBROUTINE ClumpGridInit(Info, Clump)
00148       TYPE(InfoDef) :: Info
00149       TYPE(ClumpDef), POINTER :: Clump
00150       INTEGER, POINTER, DIMENSION(:,:,:) :: mSs
00151       REAL(KIND=qPREC), POINTER, DIMENSION(:,:) :: offsets
00152       INTEGER :: nOverlaps
00153       CALL CalcPhysicalOverlaps(Info, Clump%Shape%xBounds, mSs, nOverlaps, offsets, IEVERYWHERE, lHydroPeriodic)
00154       IF (nOverlaps > 0) THEN
00155          CALL PlaceClump(Info, Clump, nOverlaps, mSs, offsets)
00156          DEALLOCATE(mSs, offsets)
00157       END IF
00158    END SUBROUTINE ClumpGridInit
00159 
00160    SUBROUTINE ClumpBeforeStep(Info, Clump)
00161       TYPE(InfoDef) :: Info
00162       TYPE(ClumpDef), POINTER :: Clump
00163       INTEGER, POINTER, DIMENSION(:,:,:) :: mSs
00164       REAL(KIND=qPREC), POINTER, DIMENSION(:,:) :: offsets
00165       INTEGER :: nOverlaps
00166       INTEGER :: i,j
00167       DO i=1,nDim
00168          DO j=1,2
00169             IF (Clump%PersistInBoundaries(i,j)) THEN
00170                CALL CalcPhysicalOverlaps(Info, Clump%Shape%xBounds, mSs, nOverlaps, offsets, IBOUNDARY(i,j), lHydroPeriodic)
00171                IF (nOverlaps > 0) THEN
00172                   CALL PlaceClump(Info, Clump, nOverlaps, mSs, offsets)
00173                   DEALLOCATE(mSs, offsets)
00174                END IF
00175             END IF
00176          END DO
00177       END DO
00178    END SUBROUTINE ClumpBeforeStep
00179 
00180    SUBROUTINE ClumpSetErrFlag(Info, Clump)
00181       TYPE(InfoDef) :: Info
00182       INTEGER, DIMENSION(3,2) :: mS    
00183       INTEGER, POINTER, DIMENSION(:,:,:) :: mSs
00184       REAL(KIND=qPREC) :: sample_res
00185       REAL(KIND=qPREC), DIMENSION(3) :: offset
00186       REAL(KIND=qPREC), POINTER, DIMENSION(:,:) :: offsets
00187       REAL(KIND=qPREC), DIMENSIOn(3,2) :: tempbounds, tempbounds2
00188       Type(ClumpDef), POINTER :: Clump
00189       INTEGER :: nOverlaps, n
00190       LOGICAL :: lrefine
00191       IF (Info%level < Clump%RefineLevel) THEN
00192          IF (levels(Info%level)%dx > .125*clump%radius) THEN !May need to trigger refinement manually            
00193             tempbounds=GetShapeBounds(Clump%Shape, levels(Info%level)%tnow)
00194             tempbounds2=GetShapeBounds(Clump%Shape, levels(Info%level)%tnow+levels(Info%level)%dt)
00195             tempbounds(:,1)=min(tempbounds(:,1), tempbounds2(:,1))
00196             tempbounds(:,2)=max(tempbounds(:,2), tempbounds2(:,2))
00197             tempbounds(1:nDim,1)=tempbounds(1:nDim,1)-.5d0*levels(Info%level)%dx
00198             tempbounds(1:nDim,2)=tempbounds(1:nDim,2)+.5d0*levels(Info%level)%dx
00199             IF (levels(Info%level)%tnow == start_time) THEN !initial setup
00200                lRefine=.true.
00201             ELSEIF((ALL(tempbounds(1:nDim,1) > GxBounds(1:nDim,1)) .AND. ALL(tempbounds(1:nDim,2) < GxBounds(1:nDim,2)))) THEN !completely inside the grid
00202                lRefine=.false.
00203             ELSEIF (ANY(tempbounds(1:nDim,2) < GxBounds(1:nDim,1)) .OR. ANY(tempbounds(1:nDim,1) > GxBounds(1:nDim,2))) THEN !completely outside the grid
00204                lRefine=.false.
00205             ELSE
00206                lRefine=.true.
00207             END IF
00208             IF (lRefine) THEN
00209                CALL CalcPhysicalOverlaps(Info, tempbounds, mSs, nOverlaps, offsets, IEVERYWHERE, lHydroPeriodic, 0)
00210                IF (nOverlaps > 0) THEN
00211                   DO n=1,nOverlaps
00212                      mS=mSs(n,:,:)
00213                      Info%ErrFlag(mS(1,1):mS(1,2), mS(2,1):mS(2,2), mS(3,1):mS(3,2))=1
00214                   END DO
00215                END IF
00216             END IF
00217          END IF
00218       END IF
00219 
00220    END SUBROUTINE ClumpSetErrFlag
00221 
00222    !includes decoding for now
00223    SUBROUTINE ClumpBeforeGlobalStep(n)
00224       INTEGER :: n
00225       !Check if clump is entirely inside of domain and if so - stop placing it in boundary zones in before step
00226       TYPE(ObjectDef), POINTER :: Object
00227       REAL(KIND=qPREC), DIMENSIOn(3,2) :: tempbounds
00228       IF (n == MaxLevel) THEN
00229          Object => ListHead
00230          DO WHILE(ASSOCIATED(Object))
00231             IF (Object%type == CLUMPOBJ) THEN
00232                pClump = transfer(Object%storage,pClump)
00233                tempbounds=GetShapeBounds(pClump%ptr%Shape, levels(MaxLevel)%tnow)
00234                IF((ALL(tempbounds(1:nDim,1) > GxBounds(1:nDim,1)) .AND. ALL(tempbounds(1:nDim,2) < GxBounds(1:nDim,2)))) THEN !completely inside the grid
00235                   pClump%ptr%PersistInBoundaries=.false.
00236                END IF
00237             END IF
00238             Object => Object%next
00239          END DO
00240       END IF
00241    END SUBROUTINE ClumpBeforeGlobalStep
00242 
00246    SUBROUTINE PlaceClump(Info,Clump, nOverlaps, mSs, offsets)
00247       TYPE(InfoDef) :: Info
00248       Type(ClumpDef) :: Clump
00249       INTEGER :: i,j,k,n,m,ii,jj,kk, location
00250       INTEGER, DIMENSION(3,2) :: mS    
00251       INTEGER, POINTER, DIMENSION(:,:,:) :: mSs
00252       REAL(KIND=qPREC), DIMENSION(3) :: offset
00253       REAL(KIND=qPREC), POINTER, DIMENSION(:,:) :: offsets
00254       REAL(KIND=qPREC), DIMENSION(3) :: xpos, pos, rpos
00255       REAL(KIND=qPREC) :: sample_fact(3), q_fact, dx,dz
00256       INTEGER :: sample_res(3), nOverlaps
00257       REAL(KIND=qPREC), DIMENSION(3,2) :: tempbounds
00258       REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: q_source, emf_source, qacc
00259       REAL(KIND=qPREC), DIMENSION(:,:,:,:), ALLOCATABLE :: emf
00260       !      RETURN
00261       dx=levels(Info%level)%dX
00262       dz=merge(dx,0d0,nDim==3)
00263       xpos=0          
00264       IF (nOverlaps > 0) THEN
00265          ALLOCATE(q_Source(NrHydroVars), qacc(NrHydroVars))
00266          sample_res=1
00267          sample_fact=0d0
00268          IF (Clump%SubSample > 0) THEN
00269             sample_res(1:nDim)=Clump%SubSample !min(Clump%SubSample,2**(MaxLevel-Info%level))
00270          ELSE
00271             sample_res(1:nDim)=1!max(Clump%SubSample,2**(MaxLevel-Info%level))
00272          END IF
00273          sample_fact(1:nDim)=1d0/REAL(sample_res(1:nDim),8)
00274          q_fact=product(sample_fact(1:nDim))
00275          DO n=1,nOverlaps
00276             mS=mSs(n,:,:)
00277             offset=offsets(n,:)
00278             !Set up emf first - aux fields first
00279             CALL ConvertTotalToInternalEnergy(Info%q(mS(1,1):mS(1,2), mS(2,1):mS(2,2),mS(3,1):mS(3,2),:))           
00280             IF (MaintainAuxArrays) THEN
00281                IF (nDim == 2) THEN
00282                   ALLOCATE(emf(mS(1,1):mS(1,2)+1, mS(2,1):mS(2,2)+1,1,3:3), emf_source(3:3))
00283                   emf=0
00284                   DO i=mS(1,1),mS(1,2)+1
00285                      xpos(1)=Info%xBounds(1,1)+offset(1)+(i-1)*dx
00286                      DO j=mS(2,1),mS(2,2)+1
00287                         xpos(2)=Info%xBounds(2,1)+offset(2)+(j-1)*dx
00288                         IF (IsInShape(Clump%Shape, xpos, rpos, levels(Info%level)%tnow)) THEN                         
00289                            emf_source(3)=emf_clump_2D(Clump, rpos) !RotateVectorFromShape(Clump%Shape, emf_2D(Clump, rpos)) !Assume shape is not rotated in 2D
00290                            emf(i,j,1,3)=emf(i,j,1,3)+emf_source(3)
00291                         END IF
00292                      END DO
00293                   END DO
00294                   CALL AddCurl(Info%aux(mS(1,1):mS(1,2)+1,mS(2,1):mS(2,2)+1,1,1:2),emf(:,:,1,3),dx)
00295                   DEALLOCATE(emf, emf_source)
00296                ELSE IF (nDim == 3) THEN
00297                   ALLOCATE(emf(mS(1,1):mS(1,2)+1, mS(2,1):mS(2,2)+1,mS(3,1):mS(3,2)+1,3), emf_source(3))
00298                   emf=0
00299                   DO i=mS(1,1),mS(1,2)
00300                      xpos(1)=Info%xBounds(1,1)+offset(1)+(i-1)*dx
00301                      DO j=mS(2,1),mS(2,2)
00302                         xpos(2)=Info%xBounds(2,1)+offset(2)+(j-1)*dx
00303                         DO k=mS(3,1),mS(3,2)
00304                            xpos(3)=Info%xBounds(3,1)+offset(3)+(k-1)*dz
00305                            DO m=1,3
00306                               pos=xpos
00307                               DO ii=1,sample_res(1)
00308                                  pos(m)=xpos(m)+(REAL(ii, 8)-half)*dx*sample_fact(m)
00309                                  IF (IsInShape(Clump%Shape, pos, rpos, levels(Info%level)%tnow)) THEN
00310                                     emf_source=RotateVectorFromShape(Clump%Shape, emf_Clump_3D(Clump, rpos))
00311                                     emf(i,j,k,m)=emf(i,j,k,m)+emf_source(m)
00312                                  END IF
00313                               END DO
00314                               emf(i,j,k,m)=emf(i,j,k,m)*sample_fact(m)
00315                            END DO
00316                         END DO
00317                      END DO
00318                      !                     CALL OutputDoubleArray(emf(i,:,:,3))
00319                   END DO
00320                   CALL AddCurl(Info%aux(mS(1,1):mS(1,2)+1,mS(2,1):mS(2,2)+1,mS(3,1):mS(3,2)+1,1:3),emf,dx)
00321                   DEALLOCATE(emf, emf_source)                 
00322                END IF
00323                CALL UpdateAux(Info, mS)          
00324             END IF
00325             ! Now set up cell centered quantities (density and momentum)
00326             DO k=mS(3,1),mS(3,2)
00327                xpos(3)=Info%xBounds(3,1)+offset(3)+REAL(k-1,8)*dz
00328                DO j=mS(2,1),mS(2,2)
00329                   xpos(2)=Info%xBounds(2,1)+offset(2)+REAL(j-1,8)*dx
00330                   DO i=mS(1,1),mS(1,2)
00331                      xpos(1)=Info%xBounds(1,1)+offset(1)+REAL(i-1,8)*dx
00332                      qacc=0
00333                      DO kk=1,sample_res(3)
00334                         pos(3)=xpos(3)+(REAL(kk, 8)-half)*dz*sample_fact(3)
00335                         DO jj=1,sample_res(2)
00336                            pos(2)=xpos(2)+(REAL(jj, 8)-half)*dx*sample_fact(2)
00337                            DO ii=1,sample_res(1)
00338                               pos(1)=xpos(1)+(REAL(ii, 8)-half)*dx*sample_fact(1)
00339                               IF (IsInShape(Clump%Shape, pos, rpos, levels(Info%level)%tnow)) THEN
00340 
00341                                  CALL q_ClumpSource(Clump,rpos,Info%q(i,j,k,1:NrHydroVars), q_source)
00342                                  qacc=qacc+q_source
00343 
00344                               END IF
00345                            END DO
00346                         END DO
00347                      END DO
00348                      Info%q(i,j,k,1:NrHydroVars)=Info%q(i,j,k,1:NrHydroVars)+qacc*q_fact                 
00349                   END DO
00350                END DO
00351             END DO
00352             CALL ConvertInternalToTotalEnergy(Info%q(mS(1,1):mS(1,2), mS(2,1):mS(2,2),mS(3,1):mS(3,2),:))           
00353          END DO
00354          DEALLOCATE(q_Source, qacc)
00355       END IF
00356 
00357    END SUBROUTINE PlaceClump
00358 
00359    SUBROUTINE q_ClumpSource(Clump, pos, q, s)
00360       TYPE(ClumpDef) :: Clump
00361       REAL(KIND=qPREC), DIMENSION(:) :: pos, q, s
00362       REAL(KIND=qPREC) :: r, f, rotation_speed(3), perturb, vel(3)
00363       r=sqrt(sum(pos(1:nDim)**2))
00364       f=smooth_tanh(r/Clump%Shape%size_param(1),Clump%thickness)
00365       perturb=Clump%Density
00366       s=0
00367       IF (Clump%density_profile==BE_profile) THEN
00368          perturb=BE_rho(REAL(r*Sqrt((4*pi*scalegrav*clump%density)/(gamma*clump%temperature))))*Clump%Density
00369       END IF
00370       IF (Clump%m2A > 0) perturb=perturb+Clump%Density*(1d0+Clump%m2A*(pos(1)**2-pos(2)**2)/r**2)  !cos(2a)=cos(a)^2-sin(a)^2
00371       IF (ASSOCIATED(Clump%DensityPerturbation)) perturb=perturb*exp(PerturbationValue(Clump%DensityPerturbation, pos))
00372 
00373       s(1)=(perturb-q(1))*f !fade into clump        
00374       vel=0
00375       IF (ASSOCIATED(Clump%VelocityPerturbation)) THEN
00376          vel=vel+VectorPerturbationValue(Clump%VelocityPerturbation,pos)
00377       END IF
00378       vel=vel*f
00379       IF (Clump%Omega /= 0) THEN 
00380          vel=vel+RotateVectorFromShape(Clump%Shape, Clump%Omega*(/pos(2),-pos(1),0d0/))
00381       END IF
00382 
00383       vel(1:nDim)=vel(1:nDim)+Clump%velocity(1:nDim)! + rotation_speed(1:nDim)
00384 
00385       s(ivx:ivx+nDim-1)=(perturb*vel(1:nDim)-q(ivx:ivx+nDim-1))*f!momentum associated with clump mass + rotational energy
00386       IF (iE .ne. 0) THEN
00387          IF (CLump%Pressure > 0) THEN 
00388             s(iE)=(gamma7*Clump%pressure-q(iE))*f !Internal energy associated with clump material
00389          ELSE
00390             s(iE)=(gamma7*perturb*Clump%temperature-q(iE))*f !Internal energy associated with clump material
00391          END IF
00392       END IF
00393       IF (Clump%iTracer .ne. 0) s(Clump%iTracer) = perturb*f-q(Clump%iTracer)
00394 
00395    END SUBROUTINE q_ClumpSource
00396 
00397 
00398    FUNCTION emf_clump_3D(Clump,pos)
00399       TYPE(ClumpDef) :: Clump
00400       REAL(KIND=qPrec) :: x_loc, s, r, w(3), emf_clump_3D(3),  tor_source, pol_source,r0,s2,rz,rz2,pos(3), f
00401       REAL(KIND=qPREC), PARAMETER :: r_m =.9 !Magnetic field cutoff in clump radii
00402       INTEGER :: i
00403       s2=pos(1)**2+pos(2)**2
00404       r=sqrt(s2+pos(3)**2)
00405       s=sqrt(s2)
00406       rz2=Clump%Shape%size_param(1)**2-pos(3)**2
00407       rz=sqrt(rz2)
00408       IF (s <= r_m*rz) THEN
00409          tor_source=Clump%B_toroidal*(r_m*rz2 - s2)/(2d0*r_m*Clump%Shape%size_param(1))
00410       ELSE
00411          tor_source=Clump%B_toroidal*(rz-s)**2/(2d0*(1d0-r_m)*Clump%Shape%size_param(1))
00412       END IF
00413       r0 = r/Clump%Shape%size_param(1)
00414       pol_source=Clump%B_poloidal/(2d0*Clump%Shape%size_param(1)**2)*(Clump%Shape%size_param(1)-r)**2  !planar version
00415       emf_clump_3D=tor_source*(/0d0,0d0,1d0/)+pol_source*(/pos(2),-pos(1),0d0/)
00416       IF (ASSOCIATED(Clump%MagneticPerturbation)) THEN         
00417          IF (r0 < r_m) THEN
00418             f=exp(2d0*r0/r_m)*(1d0-r0/r_m)**2
00419          ELSE
00420             f=0d0
00421          END IF
00422          DO i=1,3
00423             emf_clump_3D(i)=emf_clump_3D(i)+f*PerturbationValue(Clump%MagneticPerturbation(i), pos)
00424          END DO
00425       END IF
00426    END FUNCTION emf_clump_3D
00427 
00428 
00429    FUNCTION emf_clump_2D(Clump,pos)
00430       REAL(KIND=qPrec), DIMENSION(:) :: pos
00431       TYPE(ClumpDef) :: Clump
00432       REAL(KIND=qPREC) :: emf_clump_2D, r
00433       r=sqrt(SUM(pos(1:2)**2))
00434       emf_clump_2D=Clump%B_toroidal*(Clump%Shape%size_param(1)-r)
00435       IF (ASSOCIATED(Clump%MagneticPerturbation)) THEN
00436          emf_clump_2D=emf_clump_2D+PerturbationValue(Clump%MagneticPerturbation(1), pos)
00437       END IF
00438    END FUNCTION emf_clump_2D
00439 
00440 
00441 END MODULE Clumps
00442 
00443 
 All Classes Files Functions Variables