Scrambler
1
|
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