Fixed some logical errors

src/VTK/post.f90

@@ -105,7 +105,8 @@ double precision,dimension(:),allocatable::xi,yi,zi,ui,vi,wi,si,ei
 double precision,dimension(:),allocatable::s11,s12,s13,s22,s23,s33,str11,str12,str13,str22,str23,str33
 double precision,dimension(:),allocatable::azimuth1,azimuth3,dip1,dip3
 double precision zmax,dz,l1,l2,l3,n11,n12,n13,n21,n22,n23,n31,n32,n33,con,dxy
-double precision :: sselemx,sselemy,sselemz,xminls,yminls,zminls
+double precision :: sselemx,sselemy,sselemz,xminls,yminls,zminls,usum,vsum,wsum
+double precision :: uvwsum,disp,wmod
 
 double precision dxyz,x,y,z,xcut,dumpdp
 integer icut 
@@ -215,7 +216,7 @@ write(*,*) 'output v      field       ->',(output_v_field==1)
 write(*,*) 'output w      field       ->',(output_w_field==1)
 write(*,*) 'output velo   vectors     ->',(output_velo_vect==1) 
 write(*,*) 'output preiso velo vectors->',(output_preiso_velo_vect==1) 
-write(*,*) 'output iso velo vectors->',(output_preiso_velo_vect==1) 
+write(*,*) 'output iso velo vectors   ->',(output_preiso_velo_vect==1) 
 write(*,*) 'output disp   field       ->',(output_disp_field==1)
 write(*,*) 'output node press field   ->',(output_press_fieldn==1) 
 write(*,*) 'output raw press field    ->',(output_press_fielde==1) 
@@ -245,6 +246,7 @@ write(*,*) '--------------------------------------------------------------------
 
 read (7) ov%noctree,ov%nnode,ov%nleaves,ner,ov%nlsf,npcl,current_time
 nn=ov%nnode
+
 allocate(ov%on(nn))
 allocate(ov%x(nn),ov%y(nn),ov%z(nn))
 allocate(ov%u(nn),ov%v(nn),ov%w(nn),ov%wiso(nn))
@@ -280,9 +282,7 @@ read(7)(ov%x(i),ov%y(i),ov%z(i),   &
 !======================
 
 read(7) (ov%icon(1:8,ie),            &
-! Line below added by dwhipp - 12/09 
          ov%pressure(ie),            & 
-! Line below added by dwhipp - 12/09 
          ov%spressure(ie),           &
          ov%crit(ie),                &
          ov%e2d(ie),                 &
@@ -314,6 +314,7 @@ ov%nelemr=ner
 !============================
 
 read (7) (invoid(i),elvoid(i),ftr(i),rtf(i),influid(i),i=1,nn)
+
 ov%on=0
 do i=1,nn
    if (influid(i)) ov%on(i)=1
@@ -445,6 +446,11 @@ if (nest) then
   ov%x=ov%x*sselemx+xminls
   ov%y=ov%y*sselemy+yminls
   ov%z=ov%z*sselemz+zminls
+  do is=1,ov%nlsf
+    surface(is)%x=surface(is)%x*sselemx+xminls
+    surface(is)%y=surface(is)%y*sselemy+yminls
+    surface(is)%z=surface(is)%z*sselemz+zminls
+  enddo
 endif
 
 !==============================================================================
@@ -526,7 +532,7 @@ if (output_cubes==1) then
            xx=ov%x(ov%icon(k,i))
            yy=ov%y(ov%icon(k,i))
            zz=ov%z(ov%icon(k,i))
-           write(iunit,'(3e11.4)') xx,yy,zz
+           write(iunit,'(3e13.6)') xx,yy,zz
          end do
       !end if
    enddo
@@ -550,38 +556,52 @@ if (output_cubes==1) then
    write(iunit,'(a)')'SCALARS e2d float 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,ov%nleaves
-         do k=1,8 
-            write(iunit,'(e11.4)') ov%e2d(i)
+         do k=1,8
+            if (abs(ov%e2d(i)).ge.1.d99) ov%e2d(i)=sign(1.d99,ov%e2d(i))
+            if (abs(ov%e2d(i)).le.1.d-99) ov%e2d(i)=0.d0
+            write(iunit,'(e13.6)') ov%e2d(i)
          end do
    end do
    
    write(iunit,'(a)')'SCALARS w float 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,ov%nleaves
-      do k=1,8 
-      write(iunit,'(e11.4)')sum(ov%w(ov%icon(:,i)))/8.d0
+      do k=1,8
+        wsum=sum(ov%w(ov%icon(:,i)))/8.d0
+        if (abs(wsum).ge.1.d99) wsum=sign(1.d99,wsum)
+        if (abs(wsum).le.1.d-99) wsum=0.d0
+        write(iunit,'(e13.6)') wsum
       end do
    end do
    write(iunit,'(a)')'SCALARS v float 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,ov%nleaves
-      do k=1,8 
-      write(iunit,'(e11.4)')sum(ov%v(ov%icon(:,i)))/8.d0
+      do k=1,8
+        vsum=sum(ov%v(ov%icon(:,i)))/8.d0
+        if (abs(vsum).ge.1.d99) vsum=sign(1.d99,vsum)
+        if (abs(vsum).le.1.d-99) vsum=0.d0
+        write(iunit,'(e13.6)') vsum
       end do
    end do
    write(iunit,'(a)')'SCALARS u float 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,ov%nleaves
       do k=1,8 
-      write(iunit,'(e11.4)')sum(ov%u(ov%icon(:,i)))/8.d0
+        usum=sum(ov%u(ov%icon(:,i)))/8.d0
+        if (abs(usum).ge.1.d99) usum=sign(1.d99,usum)
+        if (abs(usum).le.1.d-99) usum=0.d0
+        write(iunit,'(e13.6)') usum
       end do
    end do
 
    write(iunit,'(a)')'SCALARS uvw float 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,ov%nleaves
-      do k=1,8 
-      write(iunit,'(e11.4)') sqrt((sum(ov%u(ov%icon(:,i)))/8.d0)**2+(sum(ov%v(ov%icon(:,i)))/8.d0)**2+(sum(ov%w(ov%icon(:,i)))/8.d0)**2)
+      do k=1,8
+        uvwsum=sqrt((sum(ov%u(ov%icon(:,i)))/8.d0)**2+(sum(ov%v(ov%icon(:,i)))/8.d0)**2+(sum(ov%w(ov%icon(:,i)))/8.d0)**2)
+        if (abs(uvwsum).ge.1.d99) uvwsum=sign(1.d99,uvwsum)
+        if (abs(uvwsum).le.1.d-99) uvwsum=0.d0
+        write(iunit,'(e13.6)') uvwsum
       end do
    end do
 
@@ -643,14 +663,20 @@ write(iunit,'(a11,i10)')'POINT_DATA ',nnnep
 write(iunit,'(a)')'TENSORS strain float'
 do i=1,ov%nleaves
    if (elvoid(i).eq.0) then
-   xx=sum(ov%x(ov%icon(:,i)))/8.d0
-   yy=sum(ov%y(ov%icon(:,i)))/8.d0
-   zz=sum(ov%z(ov%icon(:,i)))/8.d0
-   if (instrain(i)) then
-   write(iunit,'(3e11.4)') strain(1,:,i)
-   write(iunit,'(3e11.4)') strain(2,:,i)
-   write(iunit,'(3e11.4)') strain(3,:,i)
-   endif
+     xx=sum(ov%x(ov%icon(:,i)))/8.d0
+     yy=sum(ov%y(ov%icon(:,i)))/8.d0
+     zz=sum(ov%z(ov%icon(:,i)))/8.d0
+     if (instrain(i)) then
+       !if (strain(1,:,i).ge.1.d99) strain(1,:,i)=1.d99
+       !if (strain(1,:,i).le.1.d-99) strain(1,:,i)=1.d-99
+       !if (strain(2,:,i).ge.1.d99) strain(2,:,i)=1.d99
+       !if (strain(2,:,i).le.1.d-99) strain(2,:,i)=1.d-99
+       !if (strain(3,:,i).ge.1.d99) strain(3,:,i)=1.d99
+       !if (strain(3,:,i).le.1.d-99) strain(3,:,i)=1.d-99
+       write(iunit,'(3e13.6)') strain(1,:,i)
+       write(iunit,'(3e13.6)') strain(2,:,i)
+       write(iunit,'(3e13.6)') strain(3,:,i)
+     endif
    endif
 enddo
 close (iunit)
@@ -818,83 +844,117 @@ end do
 write(iunit,'(a11,i10)')'POINT_DATA ',nnn
 
 if(output_u_field==1) then 
-   write(iunit,'(a)')'SCALARS u float 1'
+   write(iunit,'(a)')'SCALARS u double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%u(rtf(i))
+      if (abs(ov%u(rtf(i))).ge.1.d99) ov%u(rtf(i))=sign(1.d99,ov%u(rtf(i)))
+      if (abs(ov%u(rtf(i))).le.1.d-99) ov%u(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov%u(rtf(i))
    enddo
 end if
 
 if(output_v_field==1) then 
-   write(iunit,'(a)')'SCALARS v float 1'
+   write(iunit,'(a)')'SCALARS v double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%v(rtf(i))
+      if (abs(ov%v(rtf(i))).ge.1.d99) ov%v(rtf(i))=sign(1.d99,ov%v(rtf(i)))
+      if (abs(ov%v(rtf(i))).le.1.d-99) ov%v(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov%v(rtf(i))
    enddo
 end if
 
 if(output_w_field==1) then 
-   write(iunit,'(a)')'SCALARS w float 1'
+   write(iunit,'(a)')'SCALARS w double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%w(rtf(i))
+      if (abs(ov%w(rtf(i))).ge.1.d99) ov%w(rtf(i))=sign(1.d99,ov%w(rtf(i)))
+      if (abs(ov%w(rtf(i))).le.1.d-99) ov%w(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov%w(rtf(i))
    enddo
 end if
 
-
 if(output_disp_field==1) then 
-   write(iunit,'(a)')'SCALARS disp float 1'
+   write(iunit,'(a)')'SCALARS disp double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') sqrt(ov%u(rtf(i))**2+ov%v(rtf(i))**2+ov%w(rtf(i))**2)
+      disp=sqrt(ov%u(rtf(i))**2+ov%v(rtf(i))**2+ov%w(rtf(i))**2)
+      if (abs(disp).ge.1.d99) disp=sign(1.d99,disp)
+      if (abs(disp).le.1.d-99) disp=0.d0
+      write(iunit,'(e13.6)') disp
    enddo
 end if
 
 if(output_press_fieldn==1) then 
-   write(iunit,'(a)')'SCALARS pressure_n float 1'
+   write(iunit,'(a)')'SCALARS pressure_n double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%nodal_pressure(rtf(i))
+      if (abs(ov%nodal_pressure(rtf(i))).ge.1.d99) then
+        ov%nodal_pressure(rtf(i))=sign(1.d99,ov%nodal_pressure(rtf(i)))
+      endif
+      if (abs(ov%nodal_pressure(rtf(i))).le.1.d-99) then
+        ov%nodal_pressure(rtf(i))=0.d0
+      endif
+      write(iunit,'(e13.6)') ov%nodal_pressure(rtf(i))
    enddo
 end if
 
 if(output_countnode_field==1) then 
-   write(iunit,'(a)')'SCALARS countnode float 1' 
+   write(iunit,'(a)')'SCALARS countnode double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
-   do i=1,nnn 
-      write(iunit,'(e11.4)') countnode(rtf(i)) 
+   do i=1,nnn
+      if (abs(countnode(rtf(i))).ge.1.d99) then
+        countnode(rtf(i))=sign(1.d99,countnode(rtf(i)))
+      endif
+      if (abs(countnode(rtf(i))).le.1.d-99) countnode(rtf(i))=0.d0
+      write(iunit,'(e13.6)') countnode(rtf(i))
    enddo 
 end if 
  
 if(output_e2d_fieldn==1) then 
-   write(iunit,'(a)')'SCALARS e2dn float 1'
+   write(iunit,'(a)')'SCALARS e2dn double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov_nodee2d(rtf(i))
+      if (abs(ov_nodee2d(rtf(i))).ge.1.d99) then
+        ov_nodee2d(rtf(i))=sign(1.d99,ov_nodee2d(rtf(i)))
+      endif
+      if (abs(ov_nodee2d(rtf(i))).le.1.d-99) ov_nodee2d(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov_nodee2d(rtf(i))
    enddo
 end if
 
 if(output_crit_field==1) then 
-   write(iunit,'(a)')'SCALARS crit float 1'
+   write(iunit,'(a)')'SCALARS crit double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov_nodecrit(rtf(i))
+      if (abs(ov_nodecrit(rtf(i))).ge.1.d99) then
+        ov_nodecrit(rtf(i))=sign(1.d99,ov_nodecrit(rtf(i)))
+      endif
+      if (abs(ov_nodecrit(rtf(i))).le.1.d-99) ov_nodecrit(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov_nodecrit(rtf(i))
    enddo
 end if
 
 if(output_strain_field==1) then 
-   write(iunit,'(a)')'SCALARS s float 1'
+   write(iunit,'(a)')'SCALARS s double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%strain(rtf(i))
+      if (abs(ov%strain(rtf(i))).ge.1.d99) then
+        ov%strain(rtf(i))=sign(1.d99,ov%strain(rtf(i)))
+      endif
+      if (abs(ov%strain(rtf(i))).le.1.d-99) ov%strain(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov%strain(rtf(i))
    enddo
 end if
 
 if(output_temp_field==1) then 
-   write(iunit,'(a)')'SCALARS t float 1'
+   write(iunit,'(a)')'SCALARS t double 1'
    write(iunit,'(a)')'LOOKUP_TABLE default'
    do i=1,nnn
-      write(iunit,'(e11.4)') ov%temp(rtf(i))
+      if (abs(ov%temp(rtf(i))).ge.1.d99) then
+        ov%temp(rtf(i))=sign(1.d99,ov%temp(rtf(i)))
+      endif
+      if (abs(ov%temp(rtf(i))).le.1.d-99) ov%temp(rtf(i))=0.d0
+      write(iunit,'(e13.6)') ov%temp(rtf(i))
    enddo
 end if
 
@@ -902,32 +962,51 @@ if(output_lsf_field==1) then
    do ilsf=1,ov%nlsf
       write (clsf,'(i2)') ilsf
       if (ilsf.lt.10) clsf(1:1)='0'
-      write(iunit,'(a)')'SCALARS lsf'//clsf//' float 1'
+      write(iunit,'(a)')'SCALARS lsf'//clsf//' double 1'
       write(iunit,'(a)')'LOOKUP_TABLE default'
       do i=1,nnn
-         write(iunit,'(e11.4)') ov%lsf(rtf(i),ilsf)
+         if (abs(ov%lsf(rtf(i),ilsf)).ge.1.d99) then
+           ov%lsf(rtf(i),ilsf)=sign(1.d99,ov%lsf(rtf(i),ilsf))
+         endif
+         if (abs(ov%lsf(rtf(i),ilsf)).le.1.d-99) ov%lsf(rtf(i),ilsf)=0.d0
+         write(iunit,'(e13.6)') ov%lsf(rtf(i),ilsf)
       enddo
    enddo
 end if
 
 if(output_velo_vect==1) then 
-   write(iunit,'(a)')'VECTORS velo float' 
+   write(iunit,'(a)')'VECTORS velo double' 
    do i=1,nnn 
-      write(iunit,'(3e15.4)')ov%u(rtf(i)),ov%v(rtf(i)),ov%w(rtf(i)) 
+      if (abs(ov%u(rtf(i))).ge.1.d99) ov%u(rtf(i))=sign(1.d99,ov%u(rtf(i)))
+      if (abs(ov%u(rtf(i))).le.1.d-99) ov%u(rtf(i))=0.d0
+      if (abs(ov%v(rtf(i))).ge.1.d99) ov%v(rtf(i))=sign(1.d99,ov%v(rtf(i)))
+      if (abs(ov%v(rtf(i))).le.1.d-99) ov%v(rtf(i))=0.d0
+      if (abs(ov%w(rtf(i))).ge.1.d99) ov%w(rtf(i))=sign(1.d99,ov%w(rtf(i)))
+      if (abs(ov%w(rtf(i))).le.1.d-99) ov%w(rtf(i))=0.d0
+      write(iunit,'(3e17.6)')ov%u(rtf(i)),ov%v(rtf(i)),ov%w(rtf(i))
    enddo 
 endif 
 
 if(output_preiso_velo_vect==1) then 
-   write(iunit,'(a)')'VECTORS preiso-velo float' 
+   write(iunit,'(a)')'VECTORS preiso-velo double' 
    do i=1,nnn 
-      write(iunit,'(3e15.4)')ov%u(rtf(i)),ov%v(rtf(i)),ov%w(rtf(i))-ov%wiso(rtf(i)) 
+      if (abs(ov%u(rtf(i))).ge.1.d99) ov%u(rtf(i))=sign(1.d99,ov%u(rtf(i)))
+      if (abs(ov%u(rtf(i))).le.1.d-99) ov%u(rtf(i))=0.d0
+      if (abs(ov%v(rtf(i))).ge.1.d99) ov%v(rtf(i))=sign(1.d99,ov%v(rtf(i)))
+      if (abs(ov%v(rtf(i))).le.1.d-99) ov%v(rtf(i))=0.d0
+      wmod=ov%w(rtf(i))-ov%wiso(rtf(i))
+      if (abs(wmod).ge.1.d99) wmod=sign(1.d99,wmod)
+      if (abs(wmod).le.1.d-99) wmod=0.d0
+      write(iunit,'(3e17.6)')ov%u(rtf(i)),ov%v(rtf(i)),wmod
    enddo 
 endif
 
 if(output_iso_velo_vect==1) then 
-   write(iunit,'(a)')'VECTORS iso-velo float' 
+   write(iunit,'(a)')'VECTORS iso-velo double' 
    do i=1,nnn 
-      write(iunit,'(3e15.4)')0.d0*ov%u(rtf(i)),0.d0*ov%v(rtf(i)),ov%wiso(rtf(i)) 
+      if (abs(ov%wiso(rtf(i))).ge.1.d99) ov%wiso(rtf(i))=sign(1.d99,ov%wiso(rtf(i)))
+      if (abs(ov%wiso(rtf(i))).le.1.d-99) ov%wiso(rtf(i))=0.d0
+      write(iunit,'(3e17.6)')0.d0,0.d0,ov%wiso(rtf(i))
    enddo 
 endif
 
@@ -946,11 +1025,17 @@ if (output_strain_tensor/=0) then
    do i=1,nnn
       if (nstrain(i).ne.0) strainn(:,:,i)=strainn(:,:,i)/nstrain(i)
    enddo
-   write(iunit,'(a)')'TENSORS strain float'
+   write(iunit,'(a)')'TENSORS strain double'
    do i=1,nnn
-      write(iunit,'(3e11.4)') strainn(1,:,i)
-      write(iunit,'(3e11.4)') strainn(2,:,i)
-      write(iunit,'(3e11.4)') strainn(3,:,i)
+      !if (strainn(1,:,i).ge.1.d99) strainn(1,:,i)=1.d99
+      !if (strainn(1,:,i).le.1.d-99) strainn(1,:,i)=1.d-99
+      !if (strainn(2,:,i).ge.1.d99) strainn(2,:,i)=1.d99
+      !if (strainn(2,:,i).le.1.d-99) strainn(2,:,i)=1.d-99
+      !if (strainn(3,:,i).ge.1.d99) strainn(3,:,i)=1.d99
+      !if (strainn(3,:,i).le.1.d-99) strainn(3,:,i)=1.d-99
+      write(iunit,'(3e13.6)') strainn(1,:,i)
+      write(iunit,'(3e13.6)') strainn(2,:,i)
+      write(iunit,'(3e13.6)') strainn(3,:,i)
    enddo
    deallocate (strainn,nstrain)
 endif
@@ -958,41 +1043,55 @@ endif
 write(iunit,'(a10,i10)')'CELL_DATA ', nnne 
  
 if (output_e2d_fielde==1) then
-   write(iunit,'(a)')'SCALARS e2de float 1' 
+   write(iunit,'(a)')'SCALARS e2de double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
    do ie=1,ov%nleaves 
       if (elvoid(ie).eq.0) then 
-         write(iunit,'(e11.4)') ov%e2d(ie)
+         if (abs(ov%e2d(ie)).ge.1.d99) ov%e2d(ie)=sign(1.d99,ov%e2d(ie))
+         if (abs(ov%e2d(ie)).le.1.d-99) ov%e2d(ie)=0.d0
+         write(iunit,'(e13.6)') ov%e2d(ie)
       endif 
    enddo 
 end if 
 
 if(output_press_fielde==1) then  
-   write(iunit,'(a)')'SCALARS pressure_e float 1' 
+   write(iunit,'(a)')'SCALARS pressure_e double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
    do ie=1,ov%nleaves 
       if (elvoid(ie).eq.0) then 
-         write(iunit,'(e11.4)') ov%pressure(ie) 
+         if (abs(ov%pressure(ie)).ge.1.d99) then
+           ov%pressure(ie)=sign(1.d99,ov%pressure(ie))
+         endif
+         if (abs(ov%pressure(ie)).le.1.d-99) ov%pressure(ie)=0.d0
+         write(iunit,'(e13.6)') ov%pressure(ie)
       endif 
    enddo 
 end if 
 
 if(output_smooth_press_fielde==1) then  
-   write(iunit,'(a)')'SCALARS smooth_pressure_e float 1' 
+   write(iunit,'(a)')'SCALARS smooth_pressure_e double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
    do ie=1,ov%nleaves 
       if (elvoid(ie).eq.0) then 
-         write(iunit,'(e11.4)') ov%spressure(ie)
+         if (abs(ov%spressure(ie)).ge.1.d99) then
+           ov%spressure(ie)=sign(1.d99,ov%spressure(ie))
+         endif
+         if (abs(ov%spressure(ie)).le.1.d-99) ov%spressure(ie)=0.d0
+         write(iunit,'(e13.6)') ov%spressure(ie)
       endif 
    enddo 
 end if
 
 if(output_eviscosity_fielde==1) then  
-   write(iunit,'(a)')'SCALARS eff_visc_e float 1' 
+   write(iunit,'(a)')'SCALARS eff_visc_e double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
    do ie=1,ov%nleaves 
       if (elvoid(ie).eq.0) then 
-         write(iunit,'(e11.4)') ov%eviscosity(ie) 
+         if (abs(ov%eviscosity(ie)).ge.1.d99) then
+           ov%eviscosity(ie)=sign(1.d99,ov%eviscosity(ie))
+         endif
+         if (abs(ov%eviscosity(ie)).le.1.d-99) ov%eviscosity(ie)=0.d0
+         write(iunit,'(e13.6)') ov%eviscosity(ie)
       endif 
    enddo 
 end if 
@@ -1008,11 +1107,15 @@ if(output_is_plastic_fielde==1) then
 end if
 
 if (output_dilatr_fielde==1) then
-   write(iunit,'(a)')'SCALARS dilatation_rate float 1' 
+   write(iunit,'(a)')'SCALARS dilatation_rate double 1' 
    write(iunit,'(a)')'LOOKUP_TABLE default' 
    do ie=1,ov%nleaves 
       if (elvoid(ie).eq.0) then 
-         write(iunit,'(e11.4)') ov%dilatr(ie)
+         if (abs(ov%dilatr(ie)).ge.1.d99) then
+           ov%dilatr(ie)=sign(1.d99,ov%dilatr(ie))
+         endif
+         if (abs(ov%dilatr(ie)).le.1.d-99) ov%dilatr(ie)=0.d0
+         write(iunit,'(e13.6)') ov%dilatr(ie)
       endif 
    enddo 
 end if
@@ -1061,27 +1164,27 @@ if (output_surface_icon==1) then
       write(30,'(a)')'SCALARS z float 1'
       write(30,'(a)')'LOOKUP_TABLE default'
       do i=1,surface(is)%nsurface
-         write(30,'(e11.4)') surface(is)%z(i)
+         write(30,'(e13.6)') surface(is)%z(i)
       enddo
 if (is.eq.1) then
       write(30,'(a)')'SCALARS erosion_rate float 1'
       write(30,'(a)')'LOOKUP_TABLE default'
       do i=1,surface(is)%nsurface
-         write(30,'(e11.4)') -surface(is)%r(i)
+         write(30,'(e13.6)') -surface(is)%r(i)
       enddo
       write(30,'(a)')'SCALARS rock_uplift float 1'
       write(30,'(a)')'LOOKUP_TABLE default'
       do i=1,surface(is)%nsurface
-         write(30,'(e11.4)') surface(is)%w(i)
+         write(30,'(e13.6)') surface(is)%w(i)
       enddo
 endif
       write(30,'(a)') 'VECTORS normal float'
       do i=1,surface(is)%nsurface
-         write(30,'(3e15.4)') surface(is)%xn(i),surface(is)%yn(i),surface(is)%zn(i)
+         write(30,'(3e17.6)') surface(is)%xn(i),surface(is)%yn(i),surface(is)%zn(i)
       end do
       write(30,'(a)') 'VECTORS vel float'
       do i=1,surface(is)%nsurface
-         write(30,'(3e15.4)') surface(is)%u(i),surface(is)%v(i),surface(is)%w(i)
+         write(30,'(3e17.6)') surface(is)%u(i),surface(is)%v(i),surface(is)%w(i)
       end do
       close(30)
       write(*,*) '--------------------------------------------------------------------------'
@@ -1262,38 +1365,38 @@ dip3=dip3*con
       write(30,'(a)')'SCALARS Strain float 1'
       write(30,'(a)')'LOOKUP_TABLE default'
         do i=1,(2**levmax+1)*(2**levmax+1)*nz
-        write(30,'(e11.4)') si(i)
+        write(30,'(e13.6)') si(i)
         enddo
       write(30,'(a)')'SCALARS StrainRate float 1'
       write(30,'(a)')'LOOKUP_TABLE default'
         do i=1,(2**levmax+1)*(2**levmax+1)*nz
-        write(30,'(e11.4)') ei(i)
+        write(30,'(e13.6)') ei(i)
         enddo
 !      write(30,'(a)')'SCALARS azimuth_sigma1 float 1'
 !      write(30,'(a)')'LOOKUP_TABLE default'
 !        do i=1,(2**levmax+1)*(2**levmax+1)*nz
-!        write(30,'(e11.4)') azimuth1(i)
+!        write(30,'(e13.6)') azimuth1(i)
 !        enddo
 !      write(30,'(a)')'SCALARS dip_sigma1 float 1'
 !      write(30,'(a)')'LOOKUP_TABLE default'
 !        do i=1,(2**levmax+1)*(2**levmax+1)*nz
-!        write(30,'(e11.4)') dip1(i)
+!        write(30,'(e13.6)') dip1(i)
 !        enddo
 !      write(30,'(a)')'SCALARS azimuth_sigma3 float 1'
 !      write(30,'(a)')'LOOKUP_TABLE default'
 !        do i=1,(2**levmax+1)*(2**levmax+1)*nz
-!        write(30,'(e11.4)') azimuth3(i)
+!        write(30,'(e13.6)') azimuth3(i)
 !        enddo
 !      write(30,'(a)')'SCALARS dip_sigma3 float 1'
 !      write(30,'(a)')'LOOKUP_TABLE default'
 !        do i=1,(2**levmax+1)*(2**levmax+1)*nz
-!        write(30,'(e11.4)') dip3(i)
+!        write(30,'(e13.6)') dip3(i)
 !        enddo
       write(30,'(a)')'TENSORS StrainDot float'
         do i=1,(2**levmax+1)*(2**levmax+1)*nz
-        write(30,'(3e12.4)') s11(i),s12(i),s13(i)
-        write(30,'(3e12.4)') s12(i),s22(i),s23(i)
-        write(30,'(3e12.4)') s13(i),s23(i),s33(i)
+        write(30,'(3e14.6)') s11(i),s12(i),s13(i)
+        write(30,'(3e14.6)') s12(i),s22(i),s23(i)
+        write(30,'(3e14.6)') s13(i),s23(i),s33(i)
         enddo
       close(30)
 
@@ -1342,7 +1445,7 @@ if (output_zisodisp==1) then
   write(30,'(a)')'LOOKUP_TABLE default'
   do j=1,nb+1
     do i=1,nb+1
-      write(30,'(e11.4)') zisodisp(i,j)
+      write(30,'(e13.6)') zisodisp(i,j)
     enddo
   enddo
   !write(30,'(a11,i10)') 'POINT_DATA ',int((nb+1)**2.d0)
@@ -1351,11 +1454,11 @@ if (output_zisodisp==1) then
   do j=1,nb+1
     do i=1,nb+1
       if (i==1) then
-        write(30,'(e11.4)') (zisodisp(i+1,j)-zisodisp(i,j))/dxy
+        write(30,'(e13.6)') (zisodisp(i+1,j)-zisodisp(i,j))/dxy
       elseif (i==nb+1) then
-        write(30,'(e11.4)') (zisodisp(i,j)-zisodisp(i-1,j))/dxy
+        write(30,'(e13.6)') (zisodisp(i,j)-zisodisp(i-1,j))/dxy
       else
-        write(30,'(e11.4)') ((zisodisp(i,j)-zisodisp(i-1,j))/dxy+&
+        write(30,'(e13.6)') ((zisodisp(i,j)-zisodisp(i-1,j))/dxy+&
                      (zisodisp(i+1,j)-zisodisp(i,j))/dxy)/2.d0
       endif
     enddo