From 580e34186f419a24ac72467478516d4367868ca1 Mon Sep 17 00:00:00 2001 From: davytorres Date: Mon, 20 Mar 2023 11:28:09 -0600 Subject: [PATCH] Added grid and gradient for regular communication between images --- src/matcha/grid_m.f90 | 49 +++++++++++++++++++++++++++++ src/matcha/grid_s.F90 | 70 ++++++++++++++++++++++++++++++++++++++++++ src/matcha/input_m.f90 | 44 +++++++++++++++++++++++--- src/matcha/input_s.f90 | 22 +++++++++++++ src/matcha_s.F90 | 35 ++++++++++++++++----- 5 files changed, 209 insertions(+), 11 deletions(-) create mode 100644 src/matcha/grid_m.f90 create mode 100644 src/matcha/grid_s.F90 diff --git a/src/matcha/grid_m.f90 b/src/matcha/grid_m.f90 new file mode 100644 index 0000000..d1d0bdf --- /dev/null +++ b/src/matcha/grid_m.f90 @@ -0,0 +1,49 @@ +! Copyright (c), The Regents of the University of California +! Terms of use are as specified in LICENSE.txt + +module gridr_m + implicit none + + private + public :: gridr_t + + type gridr_t + private + contains + procedure :: gridparameters + procedure :: gradient + end type + + interface gridr_t + + pure module function construct() result(gridr) + implicit none + type(gridr_t) gridr + end function + + end interface + + interface + + pure module function gridparameters(self,gb,ge,ng) result(gridp) + implicit none + class(gridr_t), intent(in) :: self + integer, intent(in) :: ng(:) + double precision, intent(in) :: gb,ge + double precision gridp(7) + end function gridparameters + + module function gradient(self, my_num_cells, ng, dx, gb, tconc, x) result(gx) + implicit none + class(gridr_t), intent(in) :: self + integer, intent(in) :: my_num_cells + integer, intent(in) :: ng(:) + double precision, intent(in) :: dx(:) + double precision, intent(in) :: gb,tconc + double precision, intent(in) :: x(:,:) + double precision, allocatable :: gx(:,:) + end function gradient + + end interface + +end module gridr_m diff --git a/src/matcha/grid_s.F90 b/src/matcha/grid_s.F90 new file mode 100644 index 0000000..762c38e --- /dev/null +++ b/src/matcha/grid_s.F90 @@ -0,0 +1,70 @@ +! Copyright (c), The Regents of the University of California +! Terms of use are as specified in LICENSE.tx +submodule(gridr_m) gridr_s + use intrinsic_array_m, only : intrinsic_array_t + +#ifdef USE_CAFFEINE + use caffeine_assert_m, only : assert +#else + use assert_m, only : assert +#endif + + implicit none + +contains + + module procedure construct + end procedure construct + + module procedure gridparameters + + gridp(1) = (ge - gb)/dble(ng(1)) + gridp(2) = (ge - gb)/dble(ng(2)) + gridp(3) = (ge - gb)/dble(ng(3)) + + gridp(4) = 1.d0/(2.d0*gridp(1)) + gridp(5) = 1.d0/(2.d0*gridp(2)) + gridp(6) = 1.d0/(2.d0*gridp(3)) + + gridp(7) = 1.d0/(gridp(1)*gridp(2)*gridp(3)) + + end procedure gridparameters + + + module procedure gradient + integer i,j,k,ii,jj,kk + double precision, allocatable :: concentration(:,:,:) + double precision, allocatable :: grad(:,:,:,:) + + allocate(concentration(ng(1),ng(2),ng(3))) + allocate(grad(ng(1),ng(2),ng(3),3)) + + concentration = 0.d0 + grad = 0.d0 + + do i = 1,my_num_cells + ii = int((x(i,1) - gb)/dx(1)) + 1 + jj = int((x(i,2) - gb)/dx(2)) + 1 + kk = int((x(i,3) - gb)/dx(3)) + 1 + concentration(ii,jj,kk) = concentration(ii,jj,kk) + tconc + end do + + call co_sum(concentration) + + do concurrent(i = 2:ng(1)-1, j = 2:ng(2)-1, k = 2:ng(3)-1) + grad(i,j,k,1) = dx(4)*(concentration(i+1,j,k)-concentration(i-1,j,k)) + grad(i,j,k,2) = dx(5)*(concentration(i,j+1,k)-concentration(i,j-1,k)) + grad(i,j,k,3) = dx(6)*(concentration(i,j,k+1)-concentration(i,j,k-1)) + end do + + allocate(gx(my_num_cells,3)) + do concurrent(i = 1:my_num_cells) + ii = int((x(i,1) - gb)/dx(1)) + 1 + jj = int((x(i,2) - gb)/dx(2)) + 1 + kk = int((x(i,3) - gb)/dx(3)) + 1 + gx(i,:) = grad(ii,jj,kk,:) + end do + + end procedure gradient + +end submodule gridr_s diff --git a/src/matcha/input_m.f90 b/src/matcha/input_m.f90 index 9d4a7c9..d0ae44c 100644 --- a/src/matcha/input_m.f90 +++ b/src/matcha/input_m.f90 @@ -8,16 +8,21 @@ module input_m type input_t private - integer :: num_cells_ = 6000, num_positions_ = 6000, num_dimensions_ = 3, num_intervals_ = 4 - double precision :: time_step_ = 0.1D0 - !double precision, allocatable :: sample_distribution_(:,:) - !allocate(sample_distribution_(num_intervals_,2)) + integer :: num_cells_ = 6000, num_positions_ = 6000, num_dimensions_ = 3, num_intervals_ = 4, ngrid_ = 50 + + double precision :: time_step_ = 0.1D-4 + double precision :: grid_begin_ = -500.d0, grid_end_ = 600.d0, cytokine_ = 1.d-6, gfac_ = 1.d-9 contains procedure :: num_cells procedure :: num_positions procedure :: num_dimensions procedure :: num_intervals + procedure :: ngrid procedure :: time_step + procedure :: grid_begin + procedure :: grid_end + procedure :: cytokine + procedure :: gfac procedure :: sample_distribution end type @@ -46,13 +51,44 @@ pure module function num_intervals(self) result(n) class(input_t), intent(in) :: self integer n end function + + pure module function ngrid(self) result(n) + implicit none + class(input_t), intent(in) :: self + integer n + end function pure module function time_step(self) result(dt) implicit none class(input_t), intent(in) :: self double precision dt end function time_step + + pure module function grid_begin(self) result(gg) + implicit none + class(input_t), intent(in) :: self + double precision gg + end function grid_begin + + pure module function grid_end(self) result(gg) + implicit none + class(input_t), intent(in) :: self + double precision gg + end function grid_end + + pure module function cytokine(self) result(gg) + implicit none + class(input_t), intent(in) :: self + double precision gg + end function cytokine + + pure module function gfac(self) result(gg) + implicit none + class(input_t), intent(in) :: self + double precision gg + end function gfac + pure module function sample_distribution(self) result(empirical_distribution) implicit none class(input_t), intent(in) :: self diff --git a/src/matcha/input_s.f90 b/src/matcha/input_s.f90 index 3492f24..a759b63 100644 --- a/src/matcha/input_s.f90 +++ b/src/matcha/input_s.f90 @@ -21,9 +21,31 @@ n = self%num_intervals_ end procedure + module procedure ngrid + n = self%ngrid_ + end procedure + module procedure time_step dt = self%time_step_ end procedure + + module procedure grid_begin + gg = self%grid_begin_ + end procedure + + module procedure grid_end + gg = self%grid_end_ + end procedure + + module procedure cytokine + gg = self%cytokine_ + end procedure + + module procedure gfac + gg = self%gfac_ + end procedure + + module procedure sample_distribution integer i,nintervals diff --git a/src/matcha_s.F90 b/src/matcha_s.F90 index 27d1c4f..d4a8c2d 100644 --- a/src/matcha_s.F90 +++ b/src/matcha_s.F90 @@ -3,8 +3,8 @@ submodule(matcha_m) matcha_s use t_cell_collection_m, only : t_cell_collection_t use distribution_m, only : distribution_t + use gridr_m, only : gridr_t use data_partition_m, only : data_partition_t - #ifdef USE_CAFFEINE use caffeine_m, only : this_image => caf_this_image #endif @@ -14,23 +14,42 @@ contains module procedure matcha - + associate( & ncells => input%num_cells(), & npositions => input%num_positions(), & ndim => input%num_dimensions(), & nintervals => input%num_intervals(), & + num_grid => input%ngrid(), & dt => input%time_step(), & + gb => input%grid_begin(), & + ge => input%grid_end(), & + cytokinev => input%cytokine(), & + gfacv => input%gfac(), & empirical_distribution => input%sample_distribution() & - ) - + ) + block double precision, parameter :: scale = 100.D0 double precision, allocatable :: random_positions(:,:), random_4vectors(:,:,:) + double precision, allocatable :: gdx(:) + double precision gxs type(distribution_t) distribution + type(gridr_t) gridr + integer, allocatable :: ng(:) integer, parameter :: nveldim = 4 + integer i,j,k,ii,jj,kk integer step + double precision tconc type(data_partition_t) data_partition + + allocate(ng(ndim)) + allocate(gdx(2*ndim+1)) + ng(:) = num_grid + gridr = gridr_t() + gdx = gridr%gridparameters(gb,ge,ng) + tconc = cytokinev*gfacv*gdx(7) + call data_partition%define_partitions(cardinality=ncells) @@ -38,7 +57,6 @@ associate(my_num_cells => data_partition%last(me) - data_partition%first(me) + 1) call random_init(repeatable=.true., image_distinct=.true.) - allocate(random_positions(my_num_cells,ndim)) call random_number(random_positions) @@ -52,8 +70,11 @@ allocate(history(nsteps)) history(1) = t_cell_collection_t(scale*random_positions, time=0.D0) do step = 2, nsteps - associate(x => history(step-1)%positions(), t => history(step-1)%time()) - history(step) = t_cell_collection_t(x + v(:,step-1,:)*dt, t + dt) + associate(x => history(step-1)%positions(), t => history(step-1)%time()) + associate(gx => gridr%gradient(my_num_cells,ng,gdx,gb,tconc,x)) + history(step) = t_cell_collection_t(x + (v(:,step-1,:) + gx(:,:))*dt, t + dt) + end associate +! history(step) = t_cell_collection_t(x + v(:,step-1,:)*dt, t + dt) end associate end do end associate