如何使用GStools生成符合指定相关函数的随机场

你使用的相关函数固定参数a=0.5时，形式为 c(x) = x^(a-2)/(x²+1)^(a/2) = x^(-1.5)/(x²+1)^0.25，可以通过GSTools的自定义协方差模型接口完成随机场生成，具体操作如下：

操作步骤 #

• 环境准备
  先安装GSTools与numpy依赖，导入所需模块：

  import numpy as np
  import gstools as gs
  

• 实现自定义协方差模型
  继承GSTools内置的CovModel基类，重写相关函数计算逻辑，同时处理x=0处的除零奇点：

  class UserDefinedCov(gs.CovModel):
      def cor(self, h):
          a = 0.5
          # 距离为0时偏移到极小值，避免数值计算报错
          h_safe = np.clip(h, 1e-12, None)
          num = np.power(h_safe, a - 2)
          den = np.power(h_safe ** 2 + 1, a / 2)
          return num / den
  

• 初始化随机场生成器
  根据你的需求指定随机场维度、方差等参数，初始化高斯随机场生成实例，固定随机种子可以复现结果：

  # 例：生成二维随机场，方差设为1，固定随机种子
  cov_model = UserDefinedCov(dim=2, var=1.0)
  srf_generator = gs.SRF(cov_model, seed=2024)
  

• 生成目标位置的随机场
  支持结构化网格、非结构化散点两种场景的随机场生成：

  # 场景1：生成100*100结构化网格上的随机场，坐标范围0~20
  x_grid = np.linspace(0, 20, 100)
  y_grid = np.linspace(0, 20, 100)
  field_structured = srf_generator((x_grid, y_grid), mesh_type="structured")
  
  # 场景2：生成1000个随机散点上的随机场
  x_scatter = np.random.uniform(0, 20, 1000)
  y_scatter = np.random.uniform(0, 20, 1000)
  field_unstructured = srf_generator((x_scatter, y_scatter))
  

注意事项 #

• 你提供的相关函数在距离x趋近于0时取值会趋向无穷大，不符合协方差函数在零距离处取有限值、随距离增大单调衰减到0的基本性质，实际运行大概率会出现数值不稳定、结果异常的问题，建议先核对相关函数的公式是否存在笔误，比如是否遗漏归一化系数、幂次位置写错等。
• 如果需要生成非高斯分布的随机场，可以在初始化SRF时传入transform参数，指定对应的分布变换规则即可。
• 如果需要调整随机场的相关长度、各向异性角度等参数，可以直接在初始化自定义协方差模型时传入对应参数，和使用GSTools内置模型的逻辑完全一致。

内容的提问来源于stack exchange，提问作者Parvin