diff --git a/pygrt/C_extension/include/grt/common/model.h b/pygrt/C_extension/include/grt/common/model.h index f22d553..ab966bb 100755 --- a/pygrt/C_extension/include/grt/common/model.h +++ b/pygrt/C_extension/include/grt/common/model.h @@ -77,12 +77,12 @@ GRT_MODEL1D * grt_init_mod1d(MYINT n); GRT_MODEL1D * grt_copy_mod1d(const GRT_MODEL1D *mod1d1); /** - * 根据不同的omega, 更新模型的 \f$ k_a^2, k_b^2 \f$ 参数 + * 根据不同的 omega, 计算衰减系数,更新弹性模量 * * @param[in,out] mod1d `MODEL1D` 结构体指针 * @param[in] omega 复数频率 */ -void grt_update_mod1d_omega(GRT_MODEL1D *mod1d, MYCOMPLEX omega); +void grt_attenuate_mod1d(GRT_MODEL1D *mod1d, MYCOMPLEX omega); /** * 根据相速度和层位,计算 iy 层的 (c/vp)^2, (c/vs)^2 以及归一化垂直波数 diff --git a/pygrt/C_extension/src/common/model.c b/pygrt/C_extension/src/common/model.c index b5b3a15..ca0c014 100755 --- a/pygrt/C_extension/src/common/model.c +++ b/pygrt/C_extension/src/common/model.c @@ -145,16 +145,17 @@ GRT_MODEL1D * grt_copy_mod1d(const GRT_MODEL1D *mod1d1){ } -void grt_update_mod1d_omega(GRT_MODEL1D *mod1d, MYCOMPLEX omega){ +void grt_attenuate_mod1d(GRT_MODEL1D *mod1d, MYCOMPLEX omega){ MYREAL Va0, Vb0; - // MYCOMPLEX ka0, kb0; MYCOMPLEX atna, atnb; for(MYINT i=0; in; ++i){ Va0 = mod1d->Va[i]; Vb0 = mod1d->Vb[i]; - atna = (mod1d->Qainv[i] > 0.0)? grt_attenuation_law(mod1d->Qainv[i], omega) : 1.0; - atnb = (mod1d->Qbinv[i] > 0.0)? grt_attenuation_law(mod1d->Qbinv[i], omega) : 1.0; + // 圆频率实部为负数表明不考虑模型的 Q 值属性 + // 在读入模型后需要需要运行一次本函数以填充弹性模量,见 grt_read_mod1d_from_file 函数 + atna = (creal(omega) >= 0.0 && mod1d->Qainv[i] > 0.0)? grt_attenuation_law(mod1d->Qainv[i], omega) : 1.0; + atnb = (creal(omega) >= 0.0 && mod1d->Qbinv[i] > 0.0)? grt_attenuation_law(mod1d->Qbinv[i], omega) : 1.0; mod1d->atna[i] = atna; mod1d->atnb[i] = atnb; @@ -391,6 +392,9 @@ GRT_MODEL1D * grt_read_mod1d_from_file(const char *command, const char *modelpat GRT_SAFE_FREE_PTR(modarr); GRT_SAFE_FREE_PTR(line); + // 填充弹性模量 + grt_attenuate_mod1d(mod1d, -1); + return mod1d; } diff --git a/pygrt/C_extension/src/dynamic/grn.c b/pygrt/C_extension/src/dynamic/grn.c index 6f4a4a3..074f329 100755 --- a/pygrt/C_extension/src/dynamic/grn.c +++ b/pygrt/C_extension/src/dynamic/grn.c @@ -159,7 +159,7 @@ void grt_integ_grn_spec( #else local_mod1d = mod1d; #endif - grt_update_mod1d_omega(local_mod1d, omega); + grt_attenuate_mod1d(local_mod1d, omega); // 是否要输出积分过程文件 bool needfstats = (statsstr!=NULL && ((grt_findElement_MYINT(statsidxs, nstatsidxs, iw) >= 0) || (grt_findElement_MYINT(statsidxs, nstatsidxs, -1) >= 0)));