femwind/fortran/ndt_boundary_conditions_test.f90

   1 program ndt_boundary_conditions_test
   2
   3 use module_boundary_conditions   ! testing only
   4 use module_utils ! to read and write matrices as text files from matlab
   5
   6 implicit none
   7
   8 real, pointer:: Kmat(:,:,:,:), &  ! fortran is not case sensitive
   9                 Kmat_m(:,:,:,:),a(:)
  10 integer :: s(4)
  11
  12 integer :: msize, &
  13     ifds, ifde, kfds, kfde, jfds, jfde,         & ! fire grid dimensions
  14     ifms, ifme, kfms, kfme, jfms, jfme,         &
  15     ifps, ifpe, kfps, kfpe, jfps, jfpe,         & ! fire patch bounds
  16     ifts, ifte, kfts, kfte, jfts, jfte
  17 integer :: i,j,k,jx
  18
  19 ! read input arrays in ikj index ordering and tight bounds
  20 call read_array_nd(a,s,'K')
  21 allocate(Kmat_m(s(1),s(2),s(3),s(4)))
  22 Kmat_m = reshape(a,s)
  23
  24 ifts = 1
  25 ifte = s(1) - 1
  26 kfts = 1
  27 kfte = s(2) - 1
  28 jfts = 1
  29 jfte = s(3) - 1
  30 ifds = 1
  31 ifde = s(1) - 1
  32 kfds = 1
  33 kfde = s(2) - 1
  34 jfds = 1
  35 jfde = s(3) - 1
  36 msize = s(4)
  37 if(msize.ne.14)call crash('msize must be 14')
  38 ifms = ifts-1
  39 ifme = ifte+2
  40 jfms = jfts-1
  41 jfme = jfte+2
  42 kfms = kfts-1
  43 kfme = kfte+2
  44
  45 ! allocate a little bigger with zeros in extra areas
  46 ! print *,'allocating Kmat size',ifms,ifme,kfms,kfme,jfms,jfme,1,msize
  47 allocate(Kmat(ifms:ifme,kfms:kfme,jfms:jfme,1:msize))
  48 Kmat = 0.
  49
  50 !print *,'copying the input matrix'
  51 Kmat(1:s(1),1:s(2),1:s(3),1:s(4))=Kmat_m(1:s(1),1:s(2),1:s(3),1:s(4))
  52
  53 !write(*,'(a)')'calling ntd_boundary_conditions'
  54 call ndt_boundary_conditions(  &
  55     ifds, ifde, kfds, kfde, jfds, jfde,         & ! fire grid dimensions
  56     ifms, ifme, kfms, kfme, jfms, jfme,         &
  57     ifps, ifpe, kfps, kfpe, jfps, jfpe,         & ! fire patch bounds
  58     ifts, ifte, kfts, kfte, jfts, jfte,         &
  59     Kmat)
  60
  61 !print *,'copying the output matrix'
  62 Kmat_m(1:s(1),1:s(2),1:s(3),1:s(4))=Kmat(1:s(1),1:s(2),1:s(3),1:s(4))
  63
  64 call write_array_nd(reshape(Kmat_m,(/product(s)/)),s,'Kb')
  65
  66 end program ndt_boundary_conditions_test