Line data Source code
1 : module coords_1d
2 :
3 : ! This module defines the Coords1D type, which is intended to to cache
4 : ! commonly used information derived from a collection of sets of 1-D
5 : ! coordinates.
6 :
7 : use shr_kind_mod, only: r8 => shr_kind_r8
8 :
9 : implicit none
10 : private
11 : save
12 :
13 : public :: Coords1D
14 :
15 : type :: Coords1D
16 : ! Number of sets of coordinates in the object.
17 : integer :: n = 0
18 : ! Number of coordinates in each set.
19 : integer :: d = 0
20 :
21 : ! All fields below will be allocated with first dimension "n".
22 : ! The second dimension is d+1 for ifc, d for mid, del, and rdel, and
23 : ! d-1 for dst and rdst.
24 :
25 : ! Cell interface coordinates.
26 : real(r8), allocatable :: ifc(:,:)
27 : ! Coordinates at cell mid-points.
28 : real(r8), allocatable :: mid(:,:)
29 : ! Width of cells.
30 : real(r8), allocatable :: del(:,:)
31 : ! Distance between cell midpoints.
32 : real(r8), allocatable :: dst(:,:)
33 : ! Reciprocals: 1/del and 1/dst.
34 : real(r8), allocatable :: rdel(:,:)
35 : real(r8), allocatable :: rdst(:,:)
36 : contains
37 : procedure :: section
38 : procedure :: finalize
39 : end type Coords1D
40 :
41 : interface Coords1D
42 : module procedure new_Coords1D_from_fields
43 : module procedure new_Coords1D_from_int
44 : end interface
45 :
46 : contains
47 :
48 : ! Constructor to create an object from existing data.
49 0 : function new_Coords1D_from_fields(ifc, mid, del, dst, &
50 0 : rdel, rdst) result(coords)
51 : real(r8), USE_CONTIGUOUS intent(in) :: ifc(:,:)
52 : real(r8), USE_CONTIGUOUS intent(in) :: mid(:,:)
53 : real(r8), USE_CONTIGUOUS intent(in) :: del(:,:)
54 : real(r8), USE_CONTIGUOUS intent(in) :: dst(:,:)
55 : real(r8), USE_CONTIGUOUS intent(in) :: rdel(:,:)
56 : real(r8), USE_CONTIGUOUS intent(in) :: rdst(:,:)
57 : type(Coords1D) :: coords
58 :
59 0 : coords = allocate_coords(size(ifc, 1), size(ifc, 2) - 1)
60 :
61 0 : coords%ifc = ifc
62 0 : coords%mid = mid
63 0 : coords%del = del
64 0 : coords%dst = dst
65 0 : coords%rdel = rdel
66 0 : coords%rdst = rdst
67 :
68 0 : end function new_Coords1D_from_fields
69 :
70 : ! Constructor if you only have interface coordinates; derives all the other
71 : ! fields.
72 8935056 : function new_Coords1D_from_int(ifc) result(coords)
73 : real(r8), USE_CONTIGUOUS intent(in) :: ifc(:,:)
74 : type(Coords1D) :: coords
75 :
76 4467528 : coords = allocate_coords(size(ifc, 1), size(ifc, 2) - 1)
77 :
78 2023062912 : coords%ifc = ifc
79 1948465584 : coords%mid = 0.5_r8 * (ifc(:,:coords%d)+ifc(:,2:))
80 1948465584 : coords%del = coords%ifc(:,2:) - coords%ifc(:,:coords%d)
81 1873868256 : coords%dst = coords%mid(:,2:) - coords%mid(:,:coords%d-1)
82 1948465584 : coords%rdel = 1._r8/coords%del
83 1873868256 : coords%rdst = 1._r8/coords%dst
84 :
85 4467528 : end function new_Coords1D_from_int
86 :
87 : ! Create a new Coords1D object that is a subsection of some other object,
88 : ! e.g. if you want only the first m coordinates, use d_bnds=[1, m].
89 : !
90 : ! Originally this used pointers, but it was found to actually be cheaper
91 : ! in practice just to make a copy, especially since pointers can impede
92 : ! optimization.
93 1489176 : function section(self, n_bnds, d_bnds)
94 : class(Coords1D), intent(in) :: self
95 : integer, intent(in) :: n_bnds(2), d_bnds(2)
96 : type(Coords1D) :: section
97 :
98 1489176 : section = allocate_coords(n_bnds(2)-n_bnds(1)+1, d_bnds(2)-d_bnds(1)+1)
99 :
100 674354304 : section%ifc = self%ifc(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2)+1)
101 649488528 : section%mid = self%mid(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2))
102 649488528 : section%del = self%del(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2))
103 624622752 : section%dst = self%dst(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2)-1)
104 649488528 : section%rdel = self%rdel(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2))
105 624622752 : section%rdst = self%rdst(n_bnds(1):n_bnds(2),d_bnds(1):d_bnds(2)-1)
106 :
107 1489176 : end function section
108 :
109 : ! Quick utility to get allocate each array with the correct size.
110 5956704 : function allocate_coords(n, d) result(coords)
111 : integer, intent(in) :: n, d
112 : type(Coords1D) :: coords
113 :
114 5956704 : coords%n = n
115 5956704 : coords%d = d
116 :
117 23826816 : allocate(coords%ifc(coords%n,coords%d+1))
118 23826816 : allocate(coords%mid(coords%n,coords%d))
119 17870112 : allocate(coords%del(coords%n,coords%d))
120 23826816 : allocate(coords%dst(coords%n,coords%d-1))
121 17870112 : allocate(coords%rdel(coords%n,coords%d))
122 17870112 : allocate(coords%rdst(coords%n,coords%d-1))
123 :
124 5956704 : end function allocate_coords
125 :
126 : ! Deallocate and reset to initial state.
127 4467528 : subroutine finalize(self)
128 : class(Coords1D), intent(inout) :: self
129 :
130 4467528 : self%n = 0
131 4467528 : self%d = 0
132 :
133 4467528 : call guarded_deallocate(self%ifc)
134 4467528 : call guarded_deallocate(self%mid)
135 4467528 : call guarded_deallocate(self%del)
136 4467528 : call guarded_deallocate(self%dst)
137 4467528 : call guarded_deallocate(self%rdel)
138 4467528 : call guarded_deallocate(self%rdst)
139 :
140 : contains
141 :
142 26805168 : subroutine guarded_deallocate(array)
143 : real(r8), allocatable :: array(:,:)
144 :
145 26805168 : if (allocated(array)) deallocate(array)
146 :
147 26805168 : end subroutine guarded_deallocate
148 :
149 : end subroutine finalize
150 :
151 0 : end module coords_1d
152 :
|