Skip to content

Curve() converts strain_unit='1' / damping_unit='1' inside the caller's array, and with interpolate=True the damping fitters get the unconverted data #85

Open

Description

@jsh9

Found while reviewing #74. Pre-existing on main (b764057); #74 does not change the behavior.

Problem

When interpolate=False, Curve.__init__() takes strain, values = data[:, 0], data[:, 1], which are views of the input array, and then converts units in place (strain *= 100; Damping_Curve also does self.damping *= 100). raw_data is the input array itself. So the caller's array is changed to %, and it is shared with the object. This is the same pattern as #59 (Ground_Motion()), here in Curve, GGmax_Curve, Damping_Curve and Stress_Curve.

When interpolate=True, the conversion happens on the interpolated copy, and raw_data stays in unit '1'. But Damping_Curve.get_HH_x_param() and get_H4_x_param() pass self.raw_data to fit_HH_x_single_layer() / fit_H4_x_single_layer(), which are documented to take both columns in %. So these fits get strain and damping values 100 times too small.

if interpolate:
strain, values = hlp.interpolate(
min_strain,
max_strain,
n_pts,
data[:, 0],
data[:, 1],
log_scale=log_scale,
)
else:
strain, values = data[:, 0], data[:, 1]
if strain_unit not in {'1', '%'}:
raise ValueError("`strain_unit` must be '1' or '%'.")
if strain_unit == '1':
strain *= 100 # strain values are internally stored in unit of %
self.raw_data = data
self.strain = strain
self.values = values

self.damping = self.values
if damping_unit not in {'1', '%'}:
raise ValueError("`damping_unit` must be '1' or '%'.")
if damping_unit == '1':
self.damping *= 100 # unit: 1 --> %

HH_x_param = hh.fit_HH_x_single_layer(
self.raw_data,
use_scipy=use_scipy,

H4_x_param = mkz.fit_H4_x_single_layer(
self.raw_data,
use_scipy=use_scipy,

damping_data_in_pct : np.ndarray
Damping data. Needs to have 2 columns (strain and damping ratio). Both
columns need to use % as unit.

damping_data_in_pct : np.ndarray
Damping data. Needs to have 2 columns (strain and damping ratio). Both
columns need to use % as unit.

from unittest import mock

import numpy as np

import PySeismoSoil.helper_hh_model as hh
from PySeismoSoil.class_curves import Damping_Curve

data_in_pct = np.genfromtxt('tests/files/curve_FKSH14.txt')[:, 2:4]  # strain [%], damping [%]

# (1) interpolate=False: the conversion to % is written into the caller's array
user_array = data_in_pct / 100.0  # strain [1], damping [1]
print(user_array[:2])
# [[1.0000e-06 1.6683e-02]
#  [3.0000e-06 1.8386e-02]]
curve = Damping_Curve(user_array, strain_unit='1', damping_unit='1')
print(user_array[:2])
# [[1.0000e-04 1.6683e+00]
#  [3.0000e-04 1.8386e+00]]  (expected unchanged)
print(curve.raw_data is user_array)
# True

# (2) interpolate=True: `raw_data` stays in unit '1', and the fitter gets it
curve_pct = Damping_Curve(data_in_pct.copy(), interpolate=True)
curve_1 = Damping_Curve(
    data_in_pct / 100.0, strain_unit='1', damping_unit='1', interpolate=True
)
for c in [curve_pct, curve_1]:
    with mock.patch.object(
        hh, 'fit_HH_x_single_layer', wraps=hh.fit_HH_x_single_layer
    ) as spy:
        c.get_HH_x_param(pop_size=1, n_gen=1)
    print(spy.call_args.args[0][:3, 0])  # strain column passed to the fitter
# [0.0001 0.0003 0.001 ]
# [1.e-06 3.e-06 1.e-05]  (expected [0.0001 0.0003 0.001 ], the same data in %)

get_H4_x_param() behaves the same: with mock.patch.object(mkz, 'fit_H4_x_single_layer', ...), the strain column passed for curve_1 is [1.e-06 3.e-06 1.e-05].

Fix

  1. In Curve.__init__(), copy the input before converting units, for example data = np.array(data, dtype=float) at the top, so the caller's array is never changed.
  2. Store raw_data in % (strain and, for Damping_Curve, damping), or have get_HH_x_param() / get_H4_x_param() build a two-column array in % from the stored data before calling the fitters. Then the fitters get % data for every strain_unit, damping_unit and interpolate combination.
  3. Add tests: the input array is unchanged after Damping_Curve(arr, strain_unit='1', damping_unit='1') and GGmax_Curve(arr, strain_unit='1'); and get_HH_x_param() / get_H4_x_param() give the same result for the same curve given in '1' and in %, with interpolate=True and interpolate=False.

Metadata

Metadata

Assignees

No one assigned

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions