Skip to content
Snippets Groups Projects
define_cloud.f90 12.3 KiB
Newer Older
  • Learn to ignore specific revisions
  • !------------------------------------------------------------------------------|
    !------------------------------------------------------------------------------|
    !                                                                              |
    !              ||===\\                                                         | 
    !              ||    \\                                                        |
    !              ||     ||   //==\\   ||  ||   //==||  ||/==\\                   |
    !              ||     ||  ||    ||  ||  ||  ||   ||  ||    ||                  |
    !              ||    //   ||    ||  ||  ||  ||   ||  ||                        |
    !              ||===//     \\==//    \\==\\  \\==\\  ||                        |
    !                                                                              |
    !------------------------------------------------------------------------------|
    !------------------------------------------------------------------------------|
    !                                                                              |
    !              DEFINE_CLOUD    Nov. 2006                                       |
    !                                                                              |
    !------------------------------------------------------------------------------|
    !------------------------------------------------------------------------------|
    
    
    subroutine define_cloud (cl,params,bcdef,nest,density_str,threadinfo)
    
    
    !------------------------------------------------------------------------------|
    !(((((((((((((((( Purpose of the routine  ))))))))))))))))))))))))))))))))))))))
    !------------------------------------------------------------------------------|
    ! if irestart=0, this routine allocates and creates the cloud of points present 
    ! in the system. Otherwise it reads from a user supplied file name the surfaces 
    ! as they were at the end of a previous run. In this case, since the run output 
    ! files contain all the octree+lsf+icloud+surface informations, the routine 
    ! first reads dummy parameters until it gets to the real interesting 
    ! cloud information 
    
    !------------------------------------------------------------------------------|
    !((((((((((((((((  declaration of the subroutine arguments  ))))))))))))))))))))
    !------------------------------------------------------------------------------|
    
    use definitions
    
    type (cloud) cl
    type (parameters) params
    
    type (nest_info) :: nest
    type (string) :: density_str(2**params%levelmax_oct,2**params%levelmax_oct)
    type (thread) threadinfo
    
    
    !------------------------------------------------------------------------------|
    !(((((((((((((((( declaration of the subroutine internal variables )))))))))))))
    !------------------------------------------------------------------------------|
    
    double precision current_time,activation_time
    double precision ssx,ssy,ssz,xls,yls,zls,eps
    logical isostasyout,isobcout,nestout,compactionout
    integer err,iproc,nproc,ierr,i,j,k,str_length
    integer ilsf,noctree,nnode,nface,nlsf,nleaves,ns,nt,restart_res
    
    integer vermajor,verminor,verstat,verrev
    
    
    !------------------------------------------------------------------------------|
    !------------------------------------------------------------------------------|
    
    call mpi_comm_size (mpi_comm_world,nproc,ierr)
    call mpi_comm_rank (mpi_comm_world,iproc,ierr)
    
    
    if (params%irestart.eq.0) then
       cl%np=0
       allocate (cl%x(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x in define_cloud$')
       allocate (cl%y(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y in define_cloud$')
       allocate (cl%z(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z in define_cloud$')
       allocate (cl%x0(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x0 in define_cloud$')
    
       allocate (cl%y0(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y0 in define_cloud$')
       allocate (cl%z0(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z0 in define_cloud$')
       allocate (cl%x0mp(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x0mp in define_cloud$')
       allocate (cl%y0mp(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y0mp in define_cloud$')
       allocate (cl%z0mp(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z0mp in define_cloud$')
    
       allocate (cl%strain(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%strain in define_cloud$')
       allocate (cl%lsf0(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%lsf0 in define_cloud$')
       allocate (cl%temp(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%temp in define_cloud$')
       allocate (cl%press(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%press in define_cloud$')
    
       allocate (cl%e2dp(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%e2dp in define_cloud$')
    
       allocate (cl%tag(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%tag in define_cloud$')
    
       allocate (cl%matnum(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%matnum in define_cloud$')
    
       allocate (cl%ematnump(1), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%ematnump in define_cloud$')
    
       if (params%isobc) then
    
         allocate (bcdef%zisodisp(2**params%levelmax_oct+1,2**params%levelmax_oct+1), stat=err); if (err.ne.0) call stop_run ('Error alloc bcdef%zisodisp in define_cloud$')
         allocate (bcdef%zisoinc(2**params%levelmax_oct+1,2**params%levelmax_oct+1), stat=err); if (err.ne.0) call stop_run ('Error alloc bcdef%zisoinc in define_cloud$')
         bcdef%zisodisp=0.d0
         bcdef%zisoinc=0.d0
    
    else
       open (19,file=trim(params%restartfile),status='old',form='unformatted')
    
       ! Read version number
       read (19) vermajor,verminor,verstat,verrev
    
    
       ! Read optional output flags
       read (19) isostasyout,isobcout,nestout,compactionout
    
    
       read (19) noctree,nnode,nleaves,nface,nlsf,cl%np,current_time     
    
    
       allocate (cl%x(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x in define_cloud$')
       allocate (cl%y(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y in define_cloud$')
       allocate (cl%z(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z in define_cloud$')
       allocate (cl%x0(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x0 in define_cloud$')
       allocate (cl%y0(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y0 in define_cloud$')
       allocate (cl%z0(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z0 in define_cloud$')
    
       allocate (cl%x0mp(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%x0mp in define_cloud$')
       allocate (cl%y0mp(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%y0mp in define_cloud$')
       allocate (cl%z0mp(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%z0mp in define_cloud$')
    
       allocate (cl%strain(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%strain in define_cloud$')
       allocate (cl%lsf0(cl%np), stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%lsf0 in define_cloud$')
       allocate (cl%temp(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%temp in define_cloud$')
       allocate (cl%press(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%press in define_cloud$')
    
       allocate (cl%e2dp(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%e2dp in define_cloud$')
    
       allocate (cl%tag(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%tag in define_cloud$')
    
       allocate (cl%matnum(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%matnum in define_cloud$')
    
       allocate (cl%ematnump(cl%np),stat=err) ; if (err.ne.0) call stop_run ('Error alloc cl%ematnump in define_cloud$')
    
         allocate (bcdef%zisodisp(2**params%levelmax_oct+1,2**params%levelmax_oct+1), stat=err); if (err.ne.0) call stop_run ('Error alloc bcdef%zisodisp in define_cloud$')
         allocate (bcdef%zisoinc(2**params%levelmax_oct+1,2**params%levelmax_oct+1), stat=err); if (err.ne.0) call stop_run ('Error alloc bcdef%zisoinc in define_cloud$')
    
       ! Dummy reads
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       if (isostasyout) read (19)
       if (compactionout) read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       if (compactionout) read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
       read (19)
    
       do ilsf=1,nlsf                                  
    
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
          read (19)
    
       read (19) cl%x
       read (19) cl%y
       read (19) cl%z
       read (19) cl%x0
       read (19) cl%y0
       read (19) cl%z0
       read (19) cl%strain
       read (19) cl%lsf0
       read (19) cl%temp
       read (19) cl%press
       read (19) cl%e2dp
       read (19) cl%tag
       read (19) cl%matnum
       read (19) cl%ematnump
    
       ! correct for vertical exaggeration
       cl%z=cl%z/params%vex
       cl%z0=cl%z0/params%vex
    
       ! read isostasy basal displacement array (optional) - dwhipp 07.13
       if (isobcout) then
         read (19) restart_res
         if (params%isobc) then
           if (restart_res /= 2**params%levelmax_oct) then
             if (iproc==0) write (*,'(a)') 'Error: Restart file and current model have different resolutions for the'
             if (iproc==0) write (*,'(a)') '       isostasy basal displacement arrays!'
             if (iproc==0) write (*,'(a)') 'This is not supported and DOUAR cannot continue.'
             call stop_run ('Error: Different octree resolutions on restart in define_cloud$')
           endif
    
         endif
       elseif (params%isobc) then
         if (iproc==0) write (*,'(a)') 'Warning: Restarting without previous isostasy basal displacement data!'
         if (iproc==0) write (*,'(a)') '         Assuming zero basal displacement.'
    
       ! read nest info (optional) - dwhipp 07.13
       if (nestout) read(19) ssx,ssy,ssz,xls,yls,zls
       if (params%nest .and. nestout) then
         if (abs(ssx-nest%sselemx) > eps .or. abs(ssy-nest%sselemy) > eps .or.     &
             abs(ssz-nest%sselemz) > eps .or. abs(xls-nest%xminls) > eps .or.      &
             abs(yls-nest%yminls) > eps .or. abs(zls-nest%zminls) > eps) then         
           if (iproc==0) write (*,'(a)') 'Error: Nest dimensions differ on restart'
           if (iproc==0) write (*,'(a)') 'This is not supported and DOUAR cannot continue.'
           call stop_run ('Error: Different nest resolutions on restart in define_cloud$')
         endif
       endif
    
       ! read compaction density strings (optional) - dwhipp 07.13
       if (compactionout) then
         read (19) restart_res
         if (params%compaction) then
           if (restart_res /= 2**params%levelmax_oct) then
             if (iproc==0) write (*,'(a)') 'Error: Restart file and current model have different resolutions for the'
             if (iproc==0) write (*,'(a)') '       compaction string arrays!'
             if (iproc==0) write (*,'(a)') 'This is not supported and DOUAR cannot continue.'
             call stop_run ('Error: Different octree resolutions on restart in define_cloud$')
           endif
           do j=1,2**params%levelmax_oct
             do i=1,2**params%levelmax_oct
               read (19) str_length
               allocate (density_str(i,j)%density(str_length),stat=threadinfo%err)
               allocate (density_str(i,j)%densityp(str_length),stat=threadinfo%err)
    
               allocate (density_str(i,j)%compactiblep(str_length),stat=threadinfo%err)
               read (19) density_str(i,j)%density,density_str(i,j)%densityp,       &
                         density_str(i,j)%compactiblep
    
             enddo
           enddo
         endif
       elseif (params%compaction) then
         if (iproc==0) write (*,'(a)') 'Warning: Restart file does not contain density string output!'
         if (iproc==0) write (*,'(a)') 'Assuming no previous sediment compaction.'
       endif
    
       close (19)
    endif
    
    if (iproc.eq.0) write (8,*) cl%np,' particles in cloud'
    
    return
    end
    
    !------------------------------------------------------------------------------|
    !------------------------------------------------------------------------------|