Added output of initial model geometry

src/DOUAR.f90

@@ -109,7 +109,7 @@ call write_streepjes(6,2)
 !call heap_init(threadinfo,2000+iproc,'./DEBUG/mpilogs/mem_heap_'//ciproc//'.dat')
 
 params%nmpi=nproc
-if (iproc.eq.0) write(*,*)  'This is DOUAR-WSMP v0.1 (2011-03-23), using ', nproc, ' processors.'
+if (iproc.eq.0) write(*,'(a,i3,a)')  'This is DOUAR-WSMP v0.1 (2011-03-23), using ', nproc, ' processors.'
 
 !------------------------------------------------------------------------------|
 !------------------------------------------------------------------------------|
@@ -125,11 +125,11 @@ numarg = 0
 
 if (numarg==0) then
    params%infile='input.txt'
-   if (iproc.eq.0) write(*,*) 'program called with no argument'
+   if (iproc.eq.0) write(*,'(a)') 'Program called with no argument'
 else
    call getarg (1,params%infile)
    if (iproc.eq.0) then
-      write(*,*) 'program called with input file ',params%infile
+      write(*,'(a,a)') 'program called with input file ',params%infile
    end if
 end if
 
@@ -744,6 +744,27 @@ do while (istep.le.params%nstep)
          ! Set reference strain rate flag to false after first nl iteration
          params%init_e2d=.false.
 
+         !--------------------------------------------------------------------------------------
+         !--------------------------------------------------------------------------------------
+         ! if this is the first iteration of the first time step, write output file
+         !--------------------------------------------------------------------------------------
+         !--------------------------------------------------------------------------------------
+
+         if (istep.eq.1 .and. iter.eq.1 .and. iter_nl.eq.1) then
+           call show_time (total,step,inc,1,'Write output initial step$')
+           do is=1,osolve%nlsf
+             allocate(surface(is)%u(surface(is)%nsurface))
+             allocate(surface(is)%v(surface(is)%nsurface))
+             allocate(surface(is)%w(surface(is)%nsurface))
+             surface(is)%u=0.d0;surface(is)%v=0.d0;surface(is)%w=0.d0
+           enddo
+           call write_global_output (params,0,0,0.d0,osolve,ov,vo,surface,cl,zi,nest,'final')
+           do is=1,osolve%nlsf
+             deallocate (surface(is)%u,surface(is)%v,surface(is)%w)
+           enddo
+           call mpi_barrier (mpi_comm_world,ierr)
+         endif
+
          !--------------------------------------------------------------------------------------
          !--------------------------------------------------------------------------------------
          ! solve system with wsmp solver