Changes between Version 6 and Version 7 of u/bliu/pnMapping

Timestamp:

11/19/21 13:37:41 (3 years ago)

Author:

Baowei Liu

Comment:

—

u/bliu/pnMapping

@@ -173 +173 @@
 
 {{{
-    ! Physical bounds of boxes
-    box_0=64000*reshape((/-1,-1,-1,1,1,1/),(/3,2/))
-    box_4=4000*reshape((/-1,-1,-1,1,1,1/),(/3,2/))
-    box_7=750*reshape((/-1,-1,-1,1,1,1/),(/3,2/))
-
-    ! index of boxes in each level's index space
-    ibox_0=nint(box_0-spread(GxBounds(:,1),2,2))/levels(0)%dx + spread((/1,0/),1,3)
-    if (maxlevel >= 4) then
-       ibox_4=nint(box_4-spread(GxBounds(:,1),2,2))/levels(4)%dx + spread((/1,0/),1,3)
-       if (maxlevel >= 7) then
-          ibox_7=nint(box_7-spread(GxBounds(:,1),2,2))/levels(7)%dx + spread((/1,0/),1,3)
-       end if
-    end if
-
-
+    nx=Info%mx(1)
+    ny=Info%mx(2)
+    nz=Info%mx(3)
+    IF (lRestart) RETURN
+    IF (.NOT. ANY(Info%level == (/0,4,7/))) RETURN
+    SELECT CASE(Info%level)
     CASE(0)
-       open(unit=195,file='NewLargeGrid_density_g_cm3_grid.dat',status='old',form='formatted',action='read')
-       open(unit=196,file='NewLargeGrid_v_x_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=197,file='NewLargeGrid_v_y_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=198,file='NewLargeGrid_v_z_cm_s_grid.dat',status='old',form='formatted',action='read')
+       open(unit=195,file='density_grid1.dat',status='old',form='formatted',action='read')
+       open(unit=196,file='v_x_grid1.dat',status='old',form='formatted',action='read')
+       open(unit=197,file='v_y_grid1.dat',status='old',form='formatted',action='read')
+       open(unit=198,file='v_z_grid1.dat',status='old',form='formatted',action='read')
        open(unit=199,file='NewLargeGrid_u_erg_g_grid.dat',status='old',form='formatted',action='read')
        ibox=ibox_0
     CASE(4)
-       open(unit=195,file='NewMedGrid_density_g_cm3_grid.dat',status='old',form='formatted',action='read')
-       open(unit=196,file='NewMedGrid_v_x_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=197,file='NewMedGrid_v_y_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=198,file='NewMedGrid_v_z_cm_s_grid.dat',status='old',form='formatted',action='read')
+       open(unit=195,file='density_grid2.dat',status='old',form='formatted',action='read')
+       open(unit=196,file='v_x_grid2.dat',status='old',form='formatted',action='read')
+       open(unit=197,file='v_y_grid2.dat',status='old',form='formatted',action='read')
+       open(unit=198,file='v_z_grid2.dat',status='old',form='formatted',action='read')
        open(unit=199,file='NewMedGrid_u_erg_g_grid.dat',status='old',form='formatted',action='read')
        ibox=ibox_4
     CASE(7)
-       open(unit=195,file='NewSmallGrid_density_g_cm3_grid.dat',status='old',form='formatted',action='read')
-       open(unit=196,file='NewSmallGrid_v_x_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=197,file='NewSmallGrid_v_y_cm_s_grid.dat',status='old',form='formatted',action='read')
-       open(unit=198,file='NewSmallGrid_v_z_cm_s_grid.dat',status='old',form='formatted',action='read')
+       open(unit=195,file='density_grid3.dat',status='old',form='formatted',action='read')
+       open(unit=196,file='v_x_grid3.dat',status='old',form='formatted',action='read')
+       open(unit=197,file='v_y_grid3.dat',status='old',form='formatted',action='read')
+       open(unit=198,file='v_z_grid3.dat',status='old',form='formatted',action='read')
        open(unit=199,file='NewSmallGrid_u_erg_g_grid.dat',status='old',form='formatted',action='read')
        ibox=ibox_7
-    select case(info%level)
-    case(0)
-        do k=1,128
-           do j=1,128
-              read(195,*),a(1,:)
-              read(196,*),a(2,:)
-              read(197,*),a(3,:)
-              read(198,*),a(4,:)
-              read(199,*),a(5,:)
-              if (j.ge.57.and.j.le.72) then
-                  do i=1,16
-                      Info%q(i,j-56,k,1)=a(1,56+i)/rscale
-                      Info%q(i,j-56,k,2)=a(1,56+i)*a(2,56+i)/rscale/velscale
-                      Info%q(i,j-56,k,3)=a(1,56+i)*a(3,56+i)/rscale/velscale
-                      Info%q(i,j-56,k,4)=a(1,56+i)*a(4,56+i)/rscale/velscale
-                      Info%q(i,j-56,k,5)=a(1,56+i)*(a(5,56+i)+0.5*(a(2,56+i)&
-                      *a(2,56+i)+a(3,56+i)*a(3,56+i)+a(4,56+i)*a(4,56+i)))/pscale
-                  end do
-              end if
-           end do
-        end do
-    case(4)
-        do k=1,128
-           do j=1,128
-              read(195,*),a(1,:)
-              read(196,*),a(2,:)
-              read(197,*),a(3,:)
-              read(198,*),a(4,:)
-              read(199,*),a(5,:)
-              Info%q(1:128,j,k,1)=a(1,:)/rscale
-              Info%q(1:128,j,k,2)=a(1,1:128)*a(2,1:128)/rscale/velscale
-              Info%q(1:128,j,k,3)=a(1,1:128)*a(3,1:128)/rscale/velscale
-              Info%q(1:128,j,k,4)=a(1,1:128)*a(4,1:128)/rscale/velscale
-              Info%q(1:128,j,k,5)=a(1,1:128)*(a(5,1:128)+0.5*(a(2,1:128)*a(2,1:128)+a(3,1:128)*a(3,1:128)+a(4,1:128)*a(4,1:128)))/pscale
-           end do
-        end do
-    case(7)
-        do k=1,192
-           do j=1,192
-              read(195,*),b(1,:)
-              read(196,*),b(2,:)
-              read(197,*),b(3,:)
-              read(198,*),b(4,:)
-              read(199,*),b(5,:)
-              Info%q(1:192,j,k,1)=b(1,:)/rscale
-              Info%q(1:192,j,k,2)=b(1,1:192)*b(2,1:192)/rscale/velscale
-              Info%q(1:192,j,k,3)=b(1,1:192)*b(3,1:192)/rscale/velscale
-              Info%q(1:192,j,k,4)=b(1,1:192)*b(4,1:192)/rscale/velscale
-              Info%q(1:192,j,k,5)=b(1,1:192)*(b(5,1:192)+0.5*(b(2,1:192)*b(2,1:192)+b(3,1:192)*b(3,1:192)+b(4,1:192)*b(4,1:192)))/pscale
-           end do
-        end do
+    END SELECT
+
+   do i=1,19  !< Reads the first 19 lines from each file - presumably containing meta data
+       read(195,*)
+       read(196,*)
+       read(197,*)
+       read(198,*)
+       read(199,*)
+    end do
+
+
+    ! Allocate space to store local box
+    allocate(a(ibox(1,1):ibox(1,2), ibox(2,1):ibox(2,2), ibox(3,1):ibox(3,2), 5))
+
+    ! Read in data
+    read(195,*),a(:,:,:,1) !< must be density in g
+    read(196,*),a(:,:,:,2) !< x-velocity in cm/s
+    read(197,*),a(:,:,:,3) !< y-velocity in cm/s
+    read(198,*),a(:,:,:,4) !< z-velocity in cm/s
+    read(199,*),a(:,:,:,5) !< internal energy in erg/g
+    close(195)
+    close(196)
+    close(197)
+    close(198)
+    close(199)
+
+    a(:,:,:,5)=a(:,:,:,5)*a(:,:,:,1)*efact / pscale !< rescale internal energy for simulation gamma
+    a(:,:,:,1)=a(:,:,:,1)/rscale
+    a(:,:,:,2:4)=a(:,:,:,2:4)/velscale * spread(a(:,:,:,1),4,3)
+    a(:,:,:,5)=a(:,:,:,5)+.5*sum(a(:,:,:,2:4)**2,4)/a(:,:,:,1)
+
+    ! Calculate overlap
+    mb(:,2)=min(ibox(:,2), info%mglobal(:,2))
+    mb(:,1)=max(ibox(:,1), info%mglobal(:,1))
+    mc=mb-spread(info%mglobal(:,1),2,2)+1
+
+    ! Update Info%q
+    Info%q(mc(1,1):mc(1,2), mc(2,1):mc(2,2), mc(3,1):mc(3,2),1:5)=a(mb(1,1):mb(1,2), mb(2,1):mb(2,2), mb(3,1):mb(3,2), :)
+    deallocate(a)
 
 }}}