Skip to content

Commit

Permalink
[WIP] Add variable pre-factors to Fourier gen.
Browse files Browse the repository at this point in the history
  • Loading branch information
LSchueler committed Mar 1, 2024
1 parent 657e8e8 commit f92fd63
Show file tree
Hide file tree
Showing 3 changed files with 25 additions and 7 deletions.
6 changes: 4 additions & 2 deletions src/gstools/field/generator.py
Original file line number Diff line number Diff line change
Expand Up @@ -635,7 +635,7 @@ def __init__(
# set model and seed
self.update(model, seed)

def __call__(self, pos, add_nugget=True):
def __call__(self, pos, add_nugget=True, phase_factor=2.*np.pi, spec_factor=1., var_factor=1.,):
"""Calculate the modes for the Fourier method.
This method calls the `summate_*` Cython methods, which are the
Expand Down Expand Up @@ -674,10 +674,12 @@ def __call__(self, pos, add_nugget=True):
self._z_1,
self._z_2,
pos,
phase_factor,
spec_factor,
)
nugget = self.get_nugget(summed_modes.shape) if add_nugget else 0.0
return (
np.sqrt(2.0 * self.model.var / np.prod(domain_size)) * summed_modes
np.sqrt(var_factor * 2.0 * self.model.var / np.prod(domain_size)) * summed_modes
+ nugget
)

Expand Down
16 changes: 15 additions & 1 deletion src/gstools/field/srf.py
Original file line number Diff line number Diff line change
Expand Up @@ -115,6 +115,9 @@ def __call__(
mesh_type="unstructured",
post_process=True,
store=True,
phase_factor=2.*np.pi,
spec_factor=1.,
var_factor=1.,
):
"""Generate the spatial random field.
Expand Down Expand Up @@ -159,7 +162,18 @@ def __call__(
# get isometrized positions and the resulting field-shape
iso_pos, shape = self.pre_pos(pos, mesh_type)
# generate the field
field = np.reshape(self.generator(iso_pos), shape)
try:
field = np.reshape(
self.generator(
iso_pos,
phase_factor=phase_factor,
spec_factor=spec_factor,
var_factor=var_factor
),
shape,
)
except TypeError:
field = np.reshape(self.generator(iso_pos, ), shape)
# upscaled variance
if not np.isscalar(point_volumes) or not np.isclose(point_volumes, 0):
scaled_var = self.upscaling_func(self.model, point_volumes)
Expand Down
10 changes: 6 additions & 4 deletions src/gstools/field/summator.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -96,7 +96,6 @@ def summate_incompr(
summed_modes[d, i] += (
proj[d] * (z_1[j] * cos(phase) + z_2[j] * sin(phase))
)

return np.asarray(summed_modes)


Expand All @@ -105,7 +104,9 @@ def summate_fourier(
const double[:, :] modes,
const double[:] z_1,
const double[:] z_2,
const double[:, :] pos
const double[:, :] pos,
const double phase_factor,
const double spec_factor,
):
cdef int i, j, d
cdef double phase
Expand All @@ -122,7 +123,8 @@ def summate_fourier(
for d in range(dim):
phase += modes[d, j] * pos[d, i]
# OpenMP doesn't like *= after +=... seems to be a compiler specific thing
phase = phase * 2. * pi
summed_modes[i] += spectral_density_sqrt[j] * (z_1[j] * cos(phase) + z_2[j] * sin(phase))
# phase = phase * 2. * pi
phase = phase * phase_factor
summed_modes[i] += spec_factor * spectral_density_sqrt[j] * (z_1[j] * cos(phase) + z_2[j] * sin(phase))

return np.asarray(summed_modes)

0 comments on commit f92fd63

Please sign in to comment.