A Modular Framework for Non-Linear Optimization of Linear Viscoelastic Model Parameters using a Second-Order Fractional Approach

<p> This repository presents a reproducible workflow for fitting the <strong>Linear Second-order Fractional (LSF)</strong> model to mastercurves of thermoplastic polymers.
The dataset is associated with the publication <em>Fauser et al. (2025)</em> and complements the dataset <a href="https://doi.org/10.18419/DARUS-4318">Data Sets for a Fractional Moisture-Dependent Viscoelasticity Model for Thermoplastic Polymers</a>.
The LSF model effectively captures the broad relaxation spectra typical of polymers, but its interdependent parameters require robust non-linear optimization. To ensure clarity and maintainability, the fitting workflow is organized into three modular components: </p> <ol> <li><strong>Viscoelastic Models Class</strong> (<code>models.py</code>)</li> <li><strong>Configuration File</strong> (<code>configuration_file.py</code>)</li> <li><strong>Model Fitting Notebook</strong> (<code>main.ipynb</code>)</li> </ol>

<h2>1. Viscoelastic Models Class (<code>models.py</code>)</h2> <p> This module serves as the core library of constitutive models. Each model is implemented as a self-contained, object-oriented class. </p> <h3>Purpose</h3> <p>Encapsulate mathematical definitions of viscoelastic models, enabling modular extensions and reuse.</p> <h3>Key Component: The LSF Model</h3> <p> The <strong>Second-Order Fractional (LSF)</strong> model is defined by its complex modulus:

</p> <pre> E*(ω) = (Ξ₀ + (iω)^α₁ Ξ₁ + (iω)^α₂ Ξ₂) / (Λ₀ + (iω)^α₁ Λ₁ + (iω)^α₂ Λ₂) </pre> <p>

Parameters:
</p> <ul> <li>Ξ₀, Ξ₁, Ξ₂ → strain-related coefficients</li>
<li>Λ₀, Λ₁, Λ₂ → stress-related coefficients</li>
<li>α₁, α₂ → fractional exponents</li> </ul>

<h3>Parent Class: <code>ViscoelasticityModel</code></h3> <p>
Defines methods for computing the real/imaginary parts of the complex modulus, as well as the storage modulus, loss modulus, and loss factor:

</p> <pre> class ViscoelasticityModel(torch.nn.Module):
"""
Parent class for viscoelasticity models
"""
def __init__(self):
super().__init__()

def complexModulus(self, w):
raise Exception

def forward(self, w):
comp = self.complexModulus(w)
return comp.imag / comp.real

def storageModulus(self,w):
comp = self.complexModulus(w)
return comp.real

def lossModulus(self,w):
comp = self.complexModulus(w)
return comp.imag</pre>

<h3>Subclass: <code>SecondOrderFractionalDerivativeModel</code></h3> <p> Implements the LSF model, with parameters represented as <code>torch.nn.Parameter</code> for optimization. A scaling factor normalizes experimental and model data for stability.

</p> <pre> class SecondOrderFractionalDerivativeModel(ViscoelasticityModel):
"""
Class for a second order derivative model
(G0 + d/dt^alpha G1 + d/dt^beta G2) stress
= (E0 + d/dt^alpha E1 + d/dt^beta E2) strain
"""
def __init__(self,E0,E1,E2,G0,G1,G2,alpha,beta):
super().__init__()

self.E0 = torch.nn.Parameter(E0 * torch.ones(1, dtype=torch.float64))
self.E1 = torch.nn.Parameter(E1 * torch.ones(1, dtype=torch.float64))
self.E2 = torch.nn.Parameter(E2 * torch.ones(1, dtype=torch.float64))

# normalize G0 to 1, G0 is not fitted
self.G0 = 1 * torch.ones(1, dtype=torch.float64)
self.G1 = torch.nn.Parameter(G1 * torch.ones(1, dtype=torch.float64))
self.G2 = torch.nn.Parameter(G2 * torch.ones(1, dtype=torch.float64))

self.alpha1 = torch.nn.Parameter(alpha * torch.ones(1, dtype=torch.float64))
self.alpha2 = torch.nn.Parameter(beta * torch.ones(1, dtype=torch.float64))

def complexModulus(self, w):

w = w[:,None]
iw = torch.view_as_complex(torch.cat((torch.zeros(w.size()),w),1))

iw_alpha1 = iw ** self.alpha1
iw_alpha2 = iw ** self.alpha2

comp = ((self.E0 + self.E1 * iw_alpha1 + self.E2 * iw_alpha2)
/ (self.G0 + self.G1 * self.E1 * iw_alpha1 + self.G2 * self.E2 * iw_alpha2))
return comp </pre>

<h2>2. Configuration File (<code>configuration_file.py</code>)</h2> <p>

This module centralizes all hyperparameters for non-linear optimization. Separating configuration ensures reproducibility, readability, and simplified parameter tuning.
</p> <h3>Example</h3>
<pre> def get_diaplex_50_optimization_params():
dict = {}
dict["lr"] = 1e-3
dict["n_iter"] = 40000
dict["data_range"] = [5e-6,1e4]
dict["w_threshold"] = 0.0032
dict["w_left"] = 1
dict["w_right"] = 2
dict["criterion"] = torch.nn.MSELoss(reduction="sum")
return dict </pre>

<h3>Key Parameters</h3> <ul> <li><code>lr</code> → learning rate for Adam optimizer</li> <li><code>n_iter</code> → number of optimization iterations</li> <li><code>data_range</code> → frequency range for fitting</li> <li><code>w_threshold</code> → frequency of peak loss factor</li> <li><code>w_left</code>, <code>w_right</code> → weighting around <code>w_threshold</code></li> <li><code>criterion</code> → loss function (L2 norm)</li> </ul>

<h2>3. Model Fitting (<code>main.ipynb</code>)</h2>

<p>
The Jupyter notebook orchestrates the full workflow for fitting <strong>fractional viscoelastic models</strong> to experimental mechanical data.
It performs data loading, preprocessing, model definition, and optimization of material parameters using <a href="https://pytorch.org/docs/stable/generated/torch.optim.Adam.html">Adam optimizer</a>.
</p>

<h3>Workflow Overview</h3>
<ol>
<li><strong>Load configuration and utilities</strong><br>
Imports helper functions and parameters from <code>configuration_file.py</code>:
<ul>
<li><code>get_diaplex_filenames()</code> — returns paths to experimental CSV data.</li>
<li><code>get_diaplex_50_optimization_params()</code> — provides hyperparameters (learning rate, iteration count, weighting, etc.).</li>
</ul>
Initializes the fractional viscoelastic model <code>SecondOrderFractionalDerivativeModel</code> for each dataset.
</li>

<li><strong>Load and preprocess experimental data</strong><br>
Reads CSV files, extracts and sorts:
<ul>
<li>Extensional Reduced Frequency</li>
<li>Extensional Storage Modulus (<em>E′</em>)</li>
<li>Extensional Loss Modulus (<em>E″</em>)</li>
</ul>
Data are converted to PyTorch tensors for GPU-based optimization.
</li>

<li><strong>Define training configuration</strong><br>
For each dataset:
<ul>
<li>Finds the frequency corresponding to the loss peak (<code>tan δ = E″/E′</code>).</li>
<li>Loads optimization settings: <code>criterion</code>, <code>lr</code>, <code>n_iter</code>, <code>w_left</code>, <code>w_right</code>, and <code>data_range</code>.</li>
<li>Computes <code>data_scaling</code> to weight frequency ranges differently around the loss peak.</li>
</ul>
</li>
</ol>

<h3>Optimization Loop</h3>
<pre>
for t in tqdm(range(n_iter)):

pred = models[key](angfreq[key][indices[key]])
meas = ((lossmod[key])/(stormod[key]))[indices[key]]
loss = criterion(meas/meas * data_scaling, pred/meas * data_scaling)

optimizer.zero_grad()
loss.backward()
optimizer.step()
</pre>

<p>
The optimization minimizes the difference between measured and predicted loss factor (<em>tan δ</em>).
Weighting (<code>w_left</code>, <code>w_right</code>) emphasizes frequencies below or above the loss peak for improved robustness.
Convergence depends strongly on initial guesses, data range, and learning rate.
</p>

<h3>Visualization and Output</h3>
<p>
After convergence, the notebook produces log–log plots comparing model predictions and measurements:
</p>
<ul>
<li><strong>Loss factor</strong> (<em>tan δ</em>)</li>
<li><strong>Storage modulus</strong> (<em>E′</em>)</li>
</ul>
<p>
Figures are saved as <code>*_fig_loss_factor.pdf</code> and <code>*_fig_storage_modulus.pdf</code>, ready for publication.
</p>

<h3>Key Notes</h3>
<ul>
<li><strong>Convergence sensitivity:</strong> Good initial parameter guesses and suitable step sizes are essential.</li>
<li><strong>Weighting strategy:</strong> Improves robustness near the loss peak, particularly for noisy data.</li>
<li><strong>Modular design:</strong> New datasets or model variants can be added easily through configuration files.</li>
</ul>

<h2>4. Fitted Parameters</h2><p>
In addition to the Python scripts for non-linear optimization, this repository also provides the fitted parameters of each master curve. These parameters, published in the related dataset
(<a href="https://doi.org/10.18419/DARUS-4318">Data Sets for a Fractional Moisture-Dependent Viscoelasticity Model for Thermoplastic Polymers</a>), are available for both the LSF model and the Generalized Maxwell model.
</p>