write_global_output.f90

!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|
!                                                                              |
!              ||===\\                                                         | 
!              ||    \\                                                        |
!              ||     ||   //==\\   ||  ||   //==||  ||/==\\                   |
!              ||     ||  ||    ||  ||  ||  ||   ||  ||    ||                  |
!              ||    //   ||    ||  ||  ||  ||   ||  ||                        |
!              ||===//     \\==//    \\==\\  \\==\\  ||                        |
!                                                                              |
!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|
!                                                                              |
!              WRITE_GLOBAL_OUTPUT     Nov. 2006                               |
!                                                                              |
!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|

subroutine write_global_output (params,istep,iter,current_time,osolve,ov, &
                                vo,surface,cl,bcdef,nest,density_str,outputtype)

!------------------------------------------------------------------------------|
!(((((((((((((((( Purpose of the routine  ))))))))))))))))))))))))))))))))))))))
!------------------------------------------------------------------------------|
! subroutine to generate large global output files in subdirectory OUT
! the files contain:
! - array lengths
! - nodal values
! - icon array
! - octree array
! - bad faces array
! - void arrays
! - surface information
! istep is time step
! osolve is solve octree
! surface are surfaces
! ns is number of surfaces
! see code for details and order of output
! Note that if you modify this routine, you may also need to modify
! define_ov.f90, define_cloud.f90, define_surface.f90 and VTK/post.f90

!------------------------------------------------------------------------------|
!((((((((((((((((  declaration of the subroutine arguments  ))))))))))))))))))))
!------------------------------------------------------------------------------|

use definitions
!use mpi

implicit none

include 'mpif.h'

type (parameters) params
integer istep,iter
double precision current_time
type (octreesolve) osolve
type (octreev) ov
type (void) vo
type (sheet) surface(params%ns)
type (cloud) cl
type (bc_definition) bcdef
type (nest_info) :: nest
type (string) :: density_str(2**params%levelmax_oct,2**params%levelmax_oct)
character*5 outputtype

!------------------------------------------------------------------------------|
!(((((((((((((((( declaration of the subroutine internal variables )))))))))))))
!------------------------------------------------------------------------------|

integer i,j,is,k,err,size_str(2),stat
integer iproc,nproc,ierr,prevstep
character*4 cistep,citer,cprevstep
character(len=5) :: outdir

!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|

call mpi_comm_size (mpi_comm_world,nproc,ierr)
call mpi_comm_rank (mpi_comm_world,iproc,ierr)

if (iproc.eq.0) then
   if (params%nest) then
     outdir=nest%ssoutdir
   else
     outdir='OUT'
   endif
   call int_to_char (cistep,4,istep)
   select case (outputtype)
      case ('debug')
         call int_to_char (citer,4,iter)
         open (9,file=trim(outdir)//'/time_'//cistep//'_'//citer//'.bin',status='unknown',form='unformatted')
      case ('final')
         open (9,file=trim(outdir)//'/time_'//cistep//'.bin',status='unknown',form='unformatted')
!         open (1234,file='OUT/time_'//cistep//'.dat',status='replace')
      case default
         call stop_run ('pb with argument in write_global_output$')
   end select

   ! write version number
   write (9) params%vermajor,params%verminor,params%verstat,params%verrev

   ! write optional output flags
   write (9) params%isostasy,(params%isobc .and. params%isostasy),params%nest, &
             params%compaction

   ! write array lengths
   write (9) osolve%octree(3),osolve%nnode,osolve%nleaves,osolve%nface,        &
             osolve%nlsf,cl%np,current_time

   ! write info on octree solve nodes (x,y,z,u,v,w,lsf,temp)
   write (9) osolve%x
   write (9) osolve%y
   write (9) osolve%z*params%vex
   write (9) osolve%u
   write (9) osolve%v
   write (9) osolve%w
   write (9) osolve%lsf
   write (9) osolve%temp
   write (9) ov%temporary_nodal_pressure
   write (9) osolve%strain
   write (9) osolve%kfix
   write (9) osolve%kfixt

   ! write isostasy and compaction output (optional) - dwhipp 07.13
   if (params%isostasy) write (9) osolve%wiso
   if (params%compaction) write (9) osolve%wcompact

   ! write icon array
   write (9) osolve%icon
   write (9) osolve%pressure
   write (9) osolve%spressure
   write (9) osolve%crit
   write (9) osolve%e2d
   write (9) osolve%eviscosity
   write (9) osolve%is_plastic
   write (9) osolve%dilatr
   write (9) osolve%matnum
   write (9) osolve%yield_ratio
   write (9) osolve%frict_angle
   write (9) ov%whole_leaf_in_fluid

   ! write compaction output (optional) - dwhipp 07.13
   if (params%compaction) write (9) osolve%compaction_density

   ! write octree information
   write (9) osolve%octree

   ! write bad face information
   write (9) osolve%iface

   ! write void information
   write (9) vo%node
   write (9) vo%leaf
   write (9) vo%ftr
   write (9) vo%rtf
   write (9) vo%influid

   ! write bad faces information
   write (9) vo%face

   ! write surface information (r,s,x,y,z,xn,yn,zn)
   do is=1,osolve%nlsf
      write (9) surface(is)%nsurface,surface(is)%activation_time,surface(is)%nt
      select case (outputtype)
      case ('debug')
        write (9) surface(is)%r
        write (9) surface(is)%s
        write (9) surface(is)%x
        write (9) surface(is)%y
        write (9) surface(is)%z*params%vex
        write (9) surface(is)%xn
        write (9) surface(is)%yn
        write (9) surface(is)%zn*params%vex
      case ('final')
        write (9) surface(is)%r
        write (9) surface(is)%s
        write (9) surface(is)%x
        write (9) surface(is)%y
        write (9) surface(is)%z*params%vex
        write (9) surface(is)%xn
        write (9) surface(is)%yn
        write (9) surface(is)%zn*params%vex
        write (9) surface(is)%u
        write (9) surface(is)%v
        write (9) surface(is)%w
      end select
      write (9) surface(is)%icon
   enddo

   ! write cloud information
   write (9) cl%x
   write (9) cl%y
   write (9) cl%z*params%vex
   write (9) cl%x0
   write (9) cl%y0
   write (9) cl%z0*params%vex
   write (9) cl%strain
   write (9) cl%lsf0
   write (9) cl%temp
   write (9) cl%press
   write (9) cl%e2dp
   write (9) cl%tag
   write (9) cl%matnum
   write (9) cl%ematnump

   ! Write isostasy basal displacement array (optional) - dwhipp 07.13
   if (params%isostasy .and. params%isobc) then
     write (9) 2**params%levelmax_oct
     write (9) bcdef%zisodisp+surface(osolve%nlsf)%sp01
!   else
!     write (9) ((0.d0+surface(osolve%nlsf)%sp01,j=1,2**params%levelmax_oct+1), &
!               i=1,2**params%levelmax_oct+1)
   end if

   ! write nested model info (optional) - dwhipp 07.13
   if (params%nest) then
     write(9) nest%sselemx,nest%sselemy,nest%sselemz,nest%xminls,nest%yminls,&
              nest%zminls
!   else
!     write(9) 0.d0,0.d0,0.d0,0.d0,0.d0,0.d0
   endif

   ! write compaction density strings - dwhipp 07.13
   if (params%compaction) then
     write (9) 2**params%levelmax_oct  ! Should test on read that this value is same as current model's val from input file
     do j=1,2**params%levelmax_oct
       do i=1,2**params%levelmax_oct
         size_str=size(density_str(i,j)%density)
         write (9) size_str(1)
         write (9) density_str(i,j)%density,density_str(i,j)%densityp,         &
                   density_str(i,j)%compactiblep
       enddo
     enddo
   endif

   close (9)

   !delete file from previous timestep if not marked for permanent save
   !mschmiddunser 12.14
   if (outputtype=='final') then
     prevstep = istep-1
       if (prevstep < params%savoffset .AND. prevstep >=0 .OR. &
         prevstep >= params%savoffset .AND. (.NOT. mod(prevstep-params%savoffset,params%savstep)==0)) then
           call int_to_char (cprevstep,4,prevstep)
           open (91,file=trim(outdir)//'/time_'//cprevstep//'.bin',iostat=stat,status='old',form='unformatted')
           if (stat.eq.0) close(91, status='delete')
       end if
   end if
end if

return
end

!------------------------------------------------------------------------------|
!------------------------------------------------------------------------------|