MODFLOW 6  version 6.6.0.dev0
USGS Modular Hydrologic Model
gwecndmodule Module Reference

Data Types

type  gwecndtype
 

Functions/Subroutines

subroutine, public cnd_cr (cndobj, name_model, input_mempath, inunit, iout, fmi, eqnsclfac, gwecommon)
 Create a new CND object. More...
 
subroutine cnd_df (this, dis, cndOptions)
 Define CND object. More...
 
subroutine cnd_ac (this, moffset, sparse)
 Add connections to CND. More...
 
subroutine cnd_mc (this, moffset, matrix_sln)
 Map CND connections. More...
 
subroutine cnd_ar (this, ibound, porosity)
 Allocate and read method for package. More...
 
subroutine cnd_ad (this)
 Advance method for the package. More...
 
subroutine cnd_fc (this, kiter, nodes, nja, matrix_sln, idxglo, rhs, cnew)
 Fill coefficient method for package. More...
 
subroutine cnd_cq (this, cnew, flowja)
 @ brief Calculate flows for package More...
 
subroutine allocate_scalars (this)
 @ brief Allocate scalar variables for package More...
 
subroutine allocate_arrays (this, nodes)
 @ brief Allocate arrays for package More...
 
subroutine cnd_da (this)
 @ brief Deallocate memory More...
 
subroutine log_options (this, found)
 Write user options to list file. More...
 
subroutine source_options (this)
 Update simulation mempath options. More...
 
subroutine log_griddata (this, found)
 Write dimensions to list file. More...
 
subroutine source_griddata (this)
 Update CND simulation data from input mempath. More...
 
subroutine calcdispellipse (this)
 Calculate dispersion coefficients. More...
 
subroutine calcdispcoef (this)
 Calculate dispersion coefficients. More...
 

Function/Subroutine Documentation

◆ allocate_arrays()

subroutine gwecndmodule::allocate_arrays ( class(gwecndtype this,
integer(i4b), intent(in)  nodes 
)

Method to allocate arrays for the package.

Definition at line 393 of file gwe-cnd.f90.

394  ! -- modules
396  use constantsmodule, only: dzero
397  ! -- dummy
398  class(GweCndType) :: this
399  integer(I4B), intent(in) :: nodes
400  !
401  ! -- Allocate
402  call mem_allocate(this%alh, nodes, 'ALH', trim(this%memoryPath))
403  call mem_allocate(this%alv, nodes, 'ALV', trim(this%memoryPath))
404  call mem_allocate(this%ath1, nodes, 'ATH1', trim(this%memoryPath))
405  call mem_allocate(this%ath2, nodes, 'ATH2', trim(this%memoryPath))
406  call mem_allocate(this%atv, nodes, 'ATV', trim(this%memoryPath))
407  call mem_allocate(this%d11, nodes, 'D11', trim(this%memoryPath))
408  call mem_allocate(this%d22, nodes, 'D22', trim(this%memoryPath))
409  call mem_allocate(this%d33, nodes, 'D33', trim(this%memoryPath))
410  call mem_allocate(this%angle1, nodes, 'ANGLE1', trim(this%memoryPath))
411  call mem_allocate(this%angle2, nodes, 'ANGLE2', trim(this%memoryPath))
412  call mem_allocate(this%angle3, nodes, 'ANGLE3', trim(this%memoryPath))
413  call mem_allocate(this%ktw, nodes, 'KTW', trim(this%memoryPath))
414  call mem_allocate(this%kts, nodes, 'KTS', trim(this%memoryPath))
415  !
416  ! -- Allocate dispersion coefficient array if xt3d not in use
417  if (this%ixt3d == 0) then
418  call mem_allocate(this%dispcoef, this%dis%njas, 'DISPCOEF', &
419  trim(this%memoryPath))
420  else
421  call mem_allocate(this%dispcoef, 0, 'DISPCOEF', trim(this%memoryPath))
422  end if
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65

◆ allocate_scalars()

subroutine gwecndmodule::allocate_scalars ( class(gwecndtype this)
private

Method to allocate scalar variables for the package.

Definition at line 342 of file gwe-cnd.f90.

343  ! -- modules
345  use constantsmodule, only: dzero
346  ! -- dummy
347  class(GweCndType) :: this
348  !
349  ! -- allocate scalars in NumericalPackageType
350  call this%NumericalPackageType%allocate_scalars()
351  !
352  ! -- Allocate
353  call mem_allocate(this%idisp, 'IDISP', this%memoryPath)
354  call mem_allocate(this%ialh, 'IALH', this%memoryPath)
355  call mem_allocate(this%ialv, 'IALV', this%memoryPath)
356  call mem_allocate(this%iath1, 'IATH1', this%memoryPath)
357  call mem_allocate(this%iath2, 'IATH2', this%memoryPath)
358  call mem_allocate(this%iatv, 'IATV', this%memoryPath)
359  call mem_allocate(this%ixt3doff, 'IXT3DOFF', this%memoryPath)
360  call mem_allocate(this%ixt3drhs, 'IXT3DRHS', this%memoryPath)
361  call mem_allocate(this%ixt3d, 'IXT3D', this%memoryPath)
362  call mem_allocate(this%id22, 'ID22', this%memoryPath)
363  call mem_allocate(this%id33, 'ID33', this%memoryPath)
364  call mem_allocate(this%iangle1, 'IANGLE1', this%memoryPath)
365  call mem_allocate(this%iangle2, 'IANGLE2', this%memoryPath)
366  call mem_allocate(this%iangle3, 'IANGLE3', this%memoryPath)
367  call mem_allocate(this%iktw, 'IKTW', this%memoryPath)
368  call mem_allocate(this%ikts, 'IKTS', this%memoryPath)
369  !
370  ! -- Initialize
371  this%idisp = 0
372  this%ialh = 0
373  this%ialv = 0
374  this%iath1 = 0
375  this%iath2 = 0
376  this%iatv = 0
377  this%ixt3doff = 0
378  this%ixt3drhs = 0
379  this%ixt3d = 0
380  this%id22 = 1
381  this%id33 = 1
382  this%iangle1 = 1
383  this%iangle2 = 1
384  this%iangle3 = 1
385  this%iktw = 1
386  this%ikts = 1

◆ calcdispcoef()

subroutine gwecndmodule::calcdispcoef ( class(gwecndtype this)
private

Definition at line 763 of file gwe-cnd.f90.

764  ! -- modules
765  use hgeoutilmodule, only: hyeff
767  ! -- dummy
768  class(GweCndType) :: this
769  ! -- local
770  integer(I4B) :: nodes, n, m, idiag, ipos
771  real(DP) :: clnm, clmn, dn, dm
772  real(DP) :: vg1, vg2, vg3
773  integer(I4B) :: ihc, isympos
774  integer(I4B) :: iavgmeth
775  real(DP) :: satn, satm, topn, topm, botn, botm
776  real(DP) :: hwva, cond, cn, cm, denom
777  real(DP) :: anm, amn, thksatn, thksatm
778  !
779  ! -- set iavgmeth = 1 to use arithmetic averaging for effective dispersion
780  iavgmeth = 1
781  !
782  ! -- Process connections
783  nodes = size(this%d11)
784  do n = 1, nodes
785  if (this%fmi%ibdgwfsat0(n) == 0) cycle
786  idiag = this%dis%con%ia(n)
787  do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
788  !
789  ! -- Set m to connected cell
790  m = this%dis%con%ja(ipos)
791  !
792  ! -- skip for lower triangle
793  if (m < n) cycle
794  isympos = this%dis%con%jas(ipos)
795  this%dispcoef(isympos) = dzero
796  if (this%fmi%ibdgwfsat0(m) == 0) cycle
797  !
798  ! -- cell dimensions
799  hwva = this%dis%con%hwva(isympos)
800  clnm = this%dis%con%cl1(isympos)
801  clmn = this%dis%con%cl2(isympos)
802  ihc = this%dis%con%ihc(isympos)
803  topn = this%dis%top(n)
804  topm = this%dis%top(m)
805  botn = this%dis%bot(n)
806  botm = this%dis%bot(m)
807  !
808  ! -- flow model information
809  satn = this%fmi%ibdgwfsat0(n)
810  satm = this%fmi%ibdgwfsat0(m)
811  !
812  ! -- Calculate dispersion coefficient for cell n in the direction
813  ! normal to the shared n-m face and for cell m in the direction
814  ! normal to the shared n-m face.
815  call this%dis%connection_normal(n, m, ihc, vg1, vg2, vg3, ipos)
816  dn = hyeff(this%d11(n), this%d22(n), this%d33(n), &
817  this%angle1(n), this%angle2(n), this%angle3(n), &
818  vg1, vg2, vg3, iavgmeth)
819  dm = hyeff(this%d11(m), this%d22(m), this%d33(m), &
820  this%angle1(m), this%angle2(m), this%angle3(m), &
821  vg1, vg2, vg3, iavgmeth)
822  !
823  ! -- Calculate dispersion conductance based on NPF subroutines and the
824  ! effective dispersion coefficients dn and dm.
825  if (ihc == 0) then
826  clnm = satn * (topn - botn) * dhalf
827  clmn = satm * (topm - botm) * dhalf
828  anm = hwva
829  !
830  ! -- n is convertible and unsaturated
831  if (satn == dzero) then
832  anm = dzero
833  else if (n > m .and. satn < done) then
834  anm = dzero
835  end if
836  !
837  ! -- m is convertible and unsaturated
838  if (satm == dzero) then
839  anm = dzero
840  else if (m > n .and. satm < done) then
841  anm = dzero
842  end if
843  !
844  ! -- amn is the same as anm for vertical flow
845  amn = anm
846  !
847  else
848  !
849  ! -- horizontal conductance
850  !
851  ! -- handle vertically staggered case
852  if (ihc == 2) then
853  thksatn = staggered_thkfrac(topn, botn, satn, topm, botm)
854  thksatm = staggered_thkfrac(topm, botm, satm, topn, botn)
855  else
856  thksatn = (topn - botn) * satn
857  thksatm = (topm - botm) * satm
858  end if
859  !
860  ! -- calculate the saturated area term
861  anm = thksatn * hwva
862  amn = thksatm * hwva
863  !
864  ! -- n or m is unsaturated, so no dispersion
865  if (satn == dzero .or. satm == dzero) then
866  anm = dzero
867  amn = dzero
868  end if
869  !
870  end if
871  !
872  ! -- calculate conductance using the two half cell conductances
873  cn = dzero
874  if (clnm > dzero) cn = dn * anm / clnm
875  cm = dzero
876  if (clmn > dzero) cm = dm * amn / clmn
877  denom = cn + cm
878  if (denom > dzero) then
879  cond = cn * cm / denom
880  else
881  cond = dzero
882  end if
883  !
884  ! -- Assign the calculated dispersion conductance
885  this%dispcoef(isympos) = cond
886  !
887  end do
888  end do
This module contains stateless conductance functions.
real(dp) function, public staggered_thkfrac(top, bot, sat, topc, botc)
Calculate the thickness fraction for staggered grids.
General-purpose hydrogeologic functions.
Definition: HGeoUtil.f90:2
real(dp) function, public hyeff(k11, k22, k33, ang1, ang2, ang3, vg1, vg2, vg3, iavgmeth)
Calculate the effective horizontal hydraulic conductivity from an ellipse using a specified direction...
Definition: HGeoUtil.f90:31
Here is the call graph for this function:

◆ calcdispellipse()

subroutine gwecndmodule::calcdispellipse ( class(gwecndtype this)

Definition at line 646 of file gwe-cnd.f90.

647  ! -- modules
648  ! -- dummy
649  class(GweCndType) :: this
650  ! -- local
651  integer(I4B) :: nodes, n
652  real(DP) :: q, qx, qy, qz
653  real(DP) :: alh, alv, ath1, ath2, atv, a
654  real(DP) :: al, at1, at2
655  real(DP) :: qzoqsquared
656  real(DP) :: ktbulk ! TODO: Implement additional options for characterizing ktbulk (see Markle refs)
657  real(DP) :: dstar
658  real(DP) :: qsw
659  !
660  ! -- loop through and calculate dispersion coefficients and angles
661  nodes = size(this%d11)
662  do n = 1, nodes
663  !
664  ! -- initialize
665  this%d11(n) = dzero
666  this%d22(n) = dzero
667  this%d33(n) = dzero
668  this%angle1(n) = dzero
669  this%angle2(n) = dzero
670  this%angle3(n) = dzero
671  if (this%fmi%ibdgwfsat0(n) == 0) cycle
672  !
673  ! -- specific discharge
674  qx = dzero
675  qy = dzero
676  qz = dzero
677  q = dzero
678  qx = this%fmi%gwfspdis(1, n)
679  qy = this%fmi%gwfspdis(2, n)
680  qz = this%fmi%gwfspdis(3, n)
681  q = qx**2 + qy**2 + qz**2
682  if (q > dzero) q = sqrt(q)
683  !
684  ! -- dispersion coefficients
685  alh = dzero
686  alv = dzero
687  ath1 = dzero
688  ath2 = dzero
689  atv = dzero
690  if (this%idisp > 0) then
691  alh = this%alh(n)
692  alv = this%alv(n)
693  ath1 = this%ath1(n)
694  ath2 = this%ath2(n)
695  atv = this%atv(n)
696  end if
697  !
698  ! -- calculate
699  ktbulk = dzero
700  if (this%iktw > 0) ktbulk = ktbulk + this%porosity(n) * this%ktw(n) * &
701  this%fmi%gwfsat(n)
702  if (this%ikts > 0) ktbulk = ktbulk + (done - this%porosity(n)) * this%kts(n)
703  !
704  ! -- The division by rhow*cpw below is undertaken to render dstar in the
705  ! form of a thermal diffusivity, and not because the governing
706  ! equation is scaled by rhow*cpw. Because of this conceptual
707  ! distinction, ktbulk is divided by the explicitly calculated product
708  ! rhow*cpw, and not by the equivalent scale factor eqnsclfac, even
709  ! though it should make no practical difference in the result.
710  dstar = ktbulk / (this%gwecommon%gwecpw * this%gwecommon%gwerhow)
711  !
712  ! -- Calculate the longitudal and transverse dispersivities
713  al = dzero
714  at1 = dzero
715  at2 = dzero
716  if (q > dzero) then
717  qzoqsquared = (qz / q)**2
718  al = alh * (done - qzoqsquared) + alv * qzoqsquared
719  at1 = ath1 * (done - qzoqsquared) + atv * qzoqsquared
720  at2 = ath2 * (done - qzoqsquared) + atv * qzoqsquared
721  end if
722  !
723  ! -- Calculate and save the diagonal components of the dispersion tensor
724  qsw = q * this%fmi%gwfsat(n) * this%eqnsclfac
725  this%d11(n) = al * qsw + ktbulk
726  this%d22(n) = at1 * qsw + ktbulk
727  this%d33(n) = at2 * qsw + ktbulk
728  !
729  ! -- Angles of rotation if velocity based dispersion tensor
730  if (this%idisp > 0) then
731  !
732  ! -- angle3 is zero
733  this%angle3(n) = dzero
734  !
735  ! -- angle2
736  a = dzero
737  if (q > dzero) a = qz / q
738  this%angle2(n) = asin(a)
739  !
740  ! -- angle1
741  a = q * cos(this%angle2(n))
742  if (a /= dzero) then
743  a = qx / a
744  else
745  a = dzero
746  end if
747  !
748  ! -- acos(1) not defined, so set to zero if necessary
749  if (a <= -done) then
750  this%angle1(n) = dpi
751  elseif (a >= done) then
752  this%angle1(n) = dzero
753  else
754  this%angle1(n) = acos(a)
755  end if
756  !
757  end if
758  end do

◆ cnd_ac()

subroutine gwecndmodule::cnd_ac ( class(gwecndtype this,
integer(i4b), intent(in)  moffset,
type(sparsematrix), intent(inout)  sparse 
)
private

Add connections for extended neighbors to the sparse matrix

Definition at line 176 of file gwe-cnd.f90.

177  ! -- modules
178  use sparsemodule, only: sparsematrix
180  ! -- dummy
181  class(GweCndType) :: this
182  integer(I4B), intent(in) :: moffset
183  type(sparsematrix), intent(inout) :: sparse
184  !
185  ! -- Add extended neighbors (neighbors of neighbors)
186  if (this%ixt3d > 0) call this%xt3d%xt3d_ac(moffset, sparse)

◆ cnd_ad()

subroutine gwecndmodule::cnd_ad ( class(gwecndtype this)
private

Definition at line 228 of file gwe-cnd.f90.

229  ! -- modules
230  use tdismodule, only: kstp, kper
231  ! -- dummy
232  class(GweCndType) :: this
233  !
234  ! -- xt3d
235  ! TODO: might consider adding a new mf6 level set pointers method, and
236  ! doing this stuff there instead of in the time step loop.
237  if (kstp * kper == 1) then
238  if (this%ixt3d > 0) then
239  call this%xt3d%xt3d_ar(this%fmi%ibdgwfsat0, this%d11, this%id33, &
240  this%d33, this%fmi%gwfsat, this%id22, this%d22, &
241  this%iangle1, this%iangle2, this%iangle3, &
242  this%angle1, this%angle2, this%angle3)
243  end if
244  end if
245  !
246  ! -- Fill d11, d22, d33, angle1, angle2, angle3 using specific discharge
247  call this%calcdispellipse()
248  !
249  ! -- Recalculate dispersion coefficients if the flows were updated
250  if (this%fmi%iflowsupdated == 1) then
251  if (this%ixt3d == 0) then
252  call this%calcdispcoef()
253  else if (this%ixt3d > 0) then
254  call this%xt3d%xt3d_fcpc(this%dis%nodes, .true.)
255  end if
256  end if
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23

◆ cnd_ar()

subroutine gwecndmodule::cnd_ar ( class(gwecndtype this,
integer(i4b), dimension(:), pointer, contiguous  ibound,
real(dp), dimension(:), pointer, contiguous  porosity 
)

Method to allocate and read static data for the package.

Definition at line 209 of file gwe-cnd.f90.

210  ! -- modules
211  ! -- dummy
212  class(GweCndType) :: this
213  integer(I4B), dimension(:), pointer, contiguous :: ibound
214  real(DP), dimension(:), pointer, contiguous :: porosity
215  ! -- local
216  ! -- formats
217  character(len=*), parameter :: fmtcnd = &
218  "(1x,/1x,'CND-- THERMAL CONDUCTION AND DISPERSION PACKAGE, VERSION 1, ', &
219  &'5/01/2023, INPUT READ FROM UNIT ', i0, //)"
220  !
221  ! -- cnd pointers to arguments that were passed in
222  this%ibound => ibound
223  this%porosity => porosity

◆ cnd_cq()

subroutine gwecndmodule::cnd_cq ( class(gwecndtype this,
real(dp), dimension(:), intent(inout)  cnew,
real(dp), dimension(:), intent(inout)  flowja 
)
private

Method to calculate dispersion contribution to flowja

Definition at line 309 of file gwe-cnd.f90.

310  ! -- modules
311  ! -- dummy
312  class(GweCndType) :: this
313  real(DP), intent(inout), dimension(:) :: cnew
314  real(DP), intent(inout), dimension(:) :: flowja
315  ! -- local
316  integer(I4B) :: n, m, ipos, isympos
317  real(DP) :: dnm, qnm
318  !
319  ! -- Calculate dispersion and add to flowja
320  if (this%ixt3d > 0) then
321  call this%xt3d%xt3d_flowja(cnew, flowja)
322  else
323  do n = 1, this%dis%nodes
324  if (this%fmi%ibdgwfsat0(n) == 0) cycle
325  do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
326  m = this%dis%con%ja(ipos)
327  if (this%fmi%ibdgwfsat0(m) == 0) cycle
328  isympos = this%dis%con%jas(ipos)
329  dnm = this%dispcoef(isympos)
330 !! qnm = dnm * (cnew(m) - cnew(n)) * this%eqnsclfac
331  qnm = dnm * (cnew(m) - cnew(n))
332  flowja(ipos) = flowja(ipos) + qnm
333  end do
334  end do
335  end if

◆ cnd_cr()

subroutine, public gwecndmodule::cnd_cr ( type(gwecndtype), pointer  cndobj,
character(len=*), intent(in)  name_model,
character(len=*), intent(in)  input_mempath,
integer(i4b), intent(in)  inunit,
integer(i4b), intent(in)  iout,
type(tspfmitype), intent(in), target  fmi,
real(dp), intent(in), pointer  eqnsclfac,
type(gweinputdatatype), intent(in), optional, target  gwecommon 
)

Create a new CND package

Parameters
[in]eqnsclfacgoverning equation scale factor
[in]gwecommonshared data container for use by multiple GWE packages

Definition at line 84 of file gwe-cnd.f90.

86  ! -- modules
87  use kindmodule, only: lgp
89  ! -- dummy
90  type(GweCndType), pointer :: cndobj
91  character(len=*), intent(in) :: name_model
92  character(len=*), intent(in) :: input_mempath
93  integer(I4B), intent(in) :: inunit
94  integer(I4B), intent(in) :: iout
95  type(TspFmiType), intent(in), target :: fmi
96  real(DP), intent(in), pointer :: eqnsclfac !< governing equation scale factor
97  type(GweInputDataType), intent(in), target, optional :: gwecommon !< shared data container for use by multiple GWE packages
98  ! -- formats
99  character(len=*), parameter :: fmtcnd = &
100  "(1x,/1x,'CND-- THERMAL CONDUCTION AND DISPERSION PACKAGE, VERSION 1, ', &
101  &'01/01/2024, INPUT READ FROM MEMPATH ', A, //)"
102  !
103  ! -- Create the object
104  allocate (cndobj)
105  !
106  ! -- create name and memory path
107  call cndobj%set_names(1, name_model, 'CND', 'CND', input_mempath)
108  !
109  ! -- Allocate scalars
110  call cndobj%allocate_scalars()
111  !
112  ! -- Set variables
113  cndobj%inunit = inunit
114  cndobj%iout = iout
115  cndobj%fmi => fmi
116  cndobj%eqnsclfac => eqnsclfac
117  if (present(gwecommon)) then
118  cndobj%gwecommon => gwecommon
119  end if
120  !
121  if (cndobj%inunit > 0) then
122  !
123  ! -- Print a message identifying the dispersion package.
124  if (cndobj%iout > 0) then
125  write (cndobj%iout, fmtcnd) input_mempath
126  end if
127  end if
This module defines variable data types.
Definition: kind.f90:8
Here is the caller graph for this function:

◆ cnd_da()

subroutine gwecndmodule::cnd_da ( class(gwecndtype this)

Method to deallocate memory for the package.

Definition at line 429 of file gwe-cnd.f90.

430  ! -- modules
434  ! -- dummy
435  class(GweCndType) :: this
436  ! -- local
437  !
438  ! -- Deallocate input memory
439  call memorystore_remove(this%name_model, 'CND', idm_context)
440  !
441  ! -- deallocate arrays
442  if (this%inunit /= 0) then
443  call mem_deallocate(this%alh, 'ALH', trim(this%memoryPath))
444  call mem_deallocate(this%alv, 'ALV', trim(this%memoryPath))
445  call mem_deallocate(this%ath1, 'ATH1', trim(this%memoryPath))
446  call mem_deallocate(this%ath2, 'ATH2', trim(this%memoryPath))
447  call mem_deallocate(this%atv, 'ATV', trim(this%memoryPath))
448  call mem_deallocate(this%d11)
449  call mem_deallocate(this%d22)
450  call mem_deallocate(this%d33)
451  call mem_deallocate(this%angle1)
452  call mem_deallocate(this%angle2)
453  call mem_deallocate(this%angle3)
454  call mem_deallocate(this%ktw)
455  call mem_deallocate(this%kts)
456  call mem_deallocate(this%dispcoef)
457  if (this%ixt3d > 0) call this%xt3d%xt3d_da()
458  end if
459  !
460  ! -- deallocate objects
461  if (this%ixt3d > 0) deallocate (this%xt3d)
462  nullify (this%gwecommon)
463  !
464  ! -- deallocate scalars
465  call mem_deallocate(this%idisp)
466  call mem_deallocate(this%ialh)
467  call mem_deallocate(this%ialv)
468  call mem_deallocate(this%iath1)
469  call mem_deallocate(this%iath2)
470  call mem_deallocate(this%iatv)
471  call mem_deallocate(this%ixt3doff)
472  call mem_deallocate(this%ixt3drhs)
473  call mem_deallocate(this%ixt3d)
474  call mem_deallocate(this%id22)
475  call mem_deallocate(this%id33)
476  call mem_deallocate(this%iangle1)
477  call mem_deallocate(this%iangle2)
478  call mem_deallocate(this%iangle3)
479  call mem_deallocate(this%iktw)
480  call mem_deallocate(this%ikts)
481  !
482  ! -- deallocate variables in NumericalPackageType
483  call this%NumericalPackageType%da()
subroutine, public memorystore_remove(component, subcomponent, context)
This module contains simulation variables.
Definition: SimVariables.f90:9
character(len=linelength) idm_context
Here is the call graph for this function:

◆ cnd_df()

subroutine gwecndmodule::cnd_df ( class(gwecndtype this,
class(disbasetype), pointer  dis,
type(gwecndoptionstype), intent(in), optional  cndOptions 
)
Parameters
[in]cndoptionsthe optional CND options, used when not creating CND from file

Definition at line 132 of file gwe-cnd.f90.

133  ! -- modules
134  ! -- dummy
135  class(GweCndType) :: this
136  class(DisBaseType), pointer :: dis
137  type(GweCndOptionsType), optional, intent(in) :: cndOptions !< the optional CND options, used when not
138  !! creating CND from file
139  !
140  ! -- Store pointer to dis
141  this%dis => dis
142  !
143  !
144  ! -- set default xt3d representation to on and lhs
145  this%ixt3d = 1
146  !
147  ! -- Read dispersion options
148  if (present(cndoptions)) then
149  this%ixt3d = cndoptions%ixt3d
150  !
151  ! -- Allocate only, grid data will not be read from file
152  call this%allocate_arrays(this%dis%nodes)
153  else
154  !
155  ! -- Source options
156  call this%source_options()
157  call this%allocate_arrays(this%dis%nodes)
158  !
159  ! -- Source dispersion data
160  call this%source_griddata()
161  end if
162  !
163  ! -- xt3d create
164  if (this%ixt3d > 0) then
165  call xt3d_cr(this%xt3d, this%name_model, this%inunit, this%iout, &
166  ldispopt=.true.)
167  this%xt3d%ixt3d = this%ixt3d
168  call this%xt3d%xt3d_df(dis)
169  end if
Here is the call graph for this function:

◆ cnd_fc()

subroutine gwecndmodule::cnd_fc ( class(gwecndtype this,
integer(i4b)  kiter,
integer(i4b), intent(in)  nodes,
integer(i4b), intent(in)  nja,
class(matrixbasetype), pointer  matrix_sln,
integer(i4b), dimension(nja), intent(in)  idxglo,
real(dp), dimension(nodes), intent(inout)  rhs,
real(dp), dimension(nodes), intent(inout)  cnew 
)

Method to calculate and fill coefficients for the package.

Definition at line 263 of file gwe-cnd.f90.

264  ! -- dummy
265  class(GweCndType) :: this
266  integer(I4B) :: kiter
267  integer(I4B), intent(in) :: nodes
268  integer(I4B), intent(in) :: nja
269  class(MatrixBaseType), pointer :: matrix_sln
270  integer(I4B), intent(in), dimension(nja) :: idxglo
271  real(DP), intent(inout), dimension(nodes) :: rhs
272  real(DP), intent(inout), dimension(nodes) :: cnew
273  ! -- local
274  integer(I4B) :: n, m, idiag, idiagm, ipos, isympos, isymcon
275  real(DP) :: dnm
276  !
277  if (this%ixt3d > 0) then
278  call this%xt3d%xt3d_fc(kiter, matrix_sln, idxglo, rhs, cnew)
279  else
280  do n = 1, nodes
281  if (this%fmi%ibdgwfsat0(n) == 0) cycle
282  idiag = this%dis%con%ia(n)
283  do ipos = this%dis%con%ia(n) + 1, this%dis%con%ia(n + 1) - 1
284  if (this%dis%con%mask(ipos) == 0) cycle
285  m = this%dis%con%ja(ipos)
286  if (m < n) cycle
287  if (this%fmi%ibdgwfsat0(m) == 0) cycle
288  isympos = this%dis%con%jas(ipos)
289  dnm = this%dispcoef(isympos)
290  !
291  ! -- Contribution to row n
292  call matrix_sln%add_value_pos(idxglo(ipos), dnm)
293  call matrix_sln%add_value_pos(idxglo(idiag), -dnm)
294  !
295  ! -- Contribution to row m
296  idiagm = this%dis%con%ia(m)
297  isymcon = this%dis%con%isym(ipos)
298  call matrix_sln%add_value_pos(idxglo(isymcon), dnm)
299  call matrix_sln%add_value_pos(idxglo(idiagm), -dnm)
300  end do
301  end do
302  end if

◆ cnd_mc()

subroutine gwecndmodule::cnd_mc ( class(gwecndtype this,
integer(i4b), intent(in)  moffset,
class(matrixbasetype), pointer  matrix_sln 
)

Map connections and construct iax, jax, and idxglox

Definition at line 193 of file gwe-cnd.f90.

194  ! -- modules
196  ! -- dummy
197  class(GweCndType) :: this
198  integer(I4B), intent(in) :: moffset
199  class(MatrixBaseType), pointer :: matrix_sln
200  !
201  ! -- Call xt3d map connections
202  if (this%ixt3d > 0) call this%xt3d%xt3d_mc(moffset, matrix_sln)

◆ log_griddata()

subroutine gwecndmodule::log_griddata ( class(gwecndtype this,
type(gwecndparamfoundtype), intent(in)  found 
)

Definition at line 528 of file gwe-cnd.f90.

530  class(GweCndType) :: this
531  type(GweCndParamFoundType), intent(in) :: found
532  !
533  write (this%iout, '(1x,a)') 'Setting CND Griddata'
534  !
535  if (found%alh) then
536  write (this%iout, '(4x,a)') 'ALH set from input file'
537  end if
538  !
539  if (found%alv) then
540  write (this%iout, '(4x,a)') 'ALV set from input file'
541  end if
542  !
543  if (found%ath1) then
544  write (this%iout, '(4x,a)') 'ATH1 set from input file'
545  end if
546  !
547  if (found%ath2) then
548  write (this%iout, '(4x,a)') 'ATH2 set from input file'
549  end if
550  !
551  if (found%atv) then
552  write (this%iout, '(4x,a)') 'ATV set from input file'
553  end if
554  !
555  if (found%ktw) then
556  write (this%iout, '(4x,a)') 'KTW set from input file'
557  end if
558  !
559  if (found%kts) then
560  write (this%iout, '(4x,a)') 'KTS set from input file'
561  end if
562  !
563  write (this%iout, '(1x,a,/)') 'End Setting CND Griddata'

◆ log_options()

subroutine gwecndmodule::log_options ( class(gwecndtype this,
type(gwecndparamfoundtype), intent(in)  found 
)

Definition at line 488 of file gwe-cnd.f90.

490  class(GweCndType) :: this
491  type(GweCndParamFoundType), intent(in) :: found
492  !
493  write (this%iout, '(1x,a)') 'Setting CND Options'
494  write (this%iout, '(4x,a,i0)') 'XT3D formulation [0=INACTIVE, 1=ACTIVE, &
495  &3=ACTIVE RHS] set to: ', this%ixt3d
496  write (this%iout, '(1x,a,/)') 'End Setting CND Options'

◆ source_griddata()

subroutine gwecndmodule::source_griddata ( class(gwecndtype this)

Definition at line 568 of file gwe-cnd.f90.

569  ! -- modules
575  ! -- dummy
576  class(GweCndType) :: this
577  ! -- locals
578  character(len=LINELENGTH) :: errmsg
579  type(GweCndParamFoundType) :: found
580  integer(I4B), dimension(:), pointer, contiguous :: map
581  ! -- formats
582  !
583  ! -- set map
584  map => null()
585  if (this%dis%nodes < this%dis%nodesuser) map => this%dis%nodeuser
586  !
587  ! -- update defaults with idm sourced values
588  call mem_set_value(this%alh, 'ALH', this%input_mempath, map, found%alh)
589  call mem_set_value(this%alv, 'ALV', this%input_mempath, map, found%alv)
590  call mem_set_value(this%ath1, 'ATH1', this%input_mempath, map, found%ath1)
591  call mem_set_value(this%ath2, 'ATH2', this%input_mempath, map, found%ath2)
592  call mem_set_value(this%atv, 'ATV', this%input_mempath, map, found%atv)
593  call mem_set_value(this%ktw, 'KTW', this%input_mempath, map, found%ktw)
594  call mem_set_value(this%kts, 'KTS', this%input_mempath, map, found%kts)
595  !
596  ! -- set active flags
597  if (found%alh) this%ialh = 1
598  if (found%alv) this%ialv = 1
599  if (found%ath1) this%iath1 = 1
600  if (found%ath2) this%iath2 = 1
601  if (found%atv) this%iatv = 1
602  if (found%ktw) this%iktw = 1
603  if (found%kts) this%ikts = 1
604  !
605  ! -- set this%idisp flag
606  if (found%alh) this%idisp = this%idisp + 1
607  if (found%alv) this%idisp = this%idisp + 1
608  if (found%ath1) this%idisp = this%idisp + 1
609  if (found%ath2) this%idisp = this%idisp + 1
610  !
611  ! -- manage dispersion arrays
612  if (this%idisp > 0) then
613  if (.not. (found%alh .and. found%ath1)) then
614  write (errmsg, '(1x,a)') &
615  'If dispersivities are specified then ALH and ATH1 are required.'
616  call store_error(errmsg)
617  end if
618  ! -- If alv not specified then point it to alh
619  if (.not. found%alv) &
620  call mem_reassignptr(this%alv, 'ALV', trim(this%memoryPath), &
621  'ALH', trim(this%memoryPath))
622  ! -- If ath2 not specified then point it to ath1
623  if (.not. found%ath2) &
624  call mem_reassignptr(this%ath2, 'ATH2', trim(this%memoryPath), &
625  'ATH1', trim(this%memoryPath))
626  ! -- If atv not specified then point it to ath2
627  if (.not. found%atv) &
628  call mem_reassignptr(this%atv, 'ATV', trim(this%memoryPath), &
629  'ATH2', trim(this%memoryPath))
630  else
631  call mem_reallocate(this%alh, 0, 'ALH', trim(this%memoryPath))
632  call mem_reallocate(this%alv, 0, 'ALV', trim(this%memoryPath))
633  call mem_reallocate(this%ath1, 0, 'ATH1', trim(this%memoryPath))
634  call mem_reallocate(this%ath2, 0, 'ATH2', trim(this%memoryPath))
635  call mem_reallocate(this%atv, 0, 'ATV', trim(this%memoryPath))
636  end if
637  !
638  ! -- log griddata
639  if (this%iout > 0) then
640  call this%log_griddata(found)
641  end if
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
integer(i4b), parameter lenmempath
maximum length of the memory path
Definition: Constants.f90:27
This module contains simulation methods.
Definition: Sim.f90:10
subroutine, public store_error(msg, terminate)
Store an error message.
Definition: Sim.f90:92
integer(i4b) function, public count_errors()
Return number of errors.
Definition: Sim.f90:59
Here is the call graph for this function:

◆ source_options()

subroutine gwecndmodule::source_options ( class(gwecndtype this)

Definition at line 501 of file gwe-cnd.f90.

502  ! -- modules
505  ! -- dummy
506  class(GweCndType) :: this
507  ! -- locals
508  type(GweCndParamFoundType) :: found
509  !
510  ! -- update defaults with idm sourced values
511  call mem_set_value(this%ixt3doff, 'XT3D_OFF', this%input_mempath, &
512  found%xt3d_off)
513  call mem_set_value(this%ixt3drhs, 'XT3D_RHS', this%input_mempath, &
514  found%xt3d_rhs)
515  !
516  ! -- set xt3d state flag
517  if (found%xt3d_off) this%ixt3d = 0
518  if (found%xt3d_rhs) this%ixt3d = 2
519  !
520  ! -- log options
521  if (this%iout > 0) then
522  call this%log_options(found)
523  end if