Skip to content

EPNN

# linux
wget -c https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstate-16-plas.dat -P ./datasets/
wget -c https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstress-16-plas.dat -P ./datasets/
# windows
# curl https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstate-16-plas.dat --create-dirs -o ./datasets/dstate-16-plas.dat
# curl https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstress-16-plas.dat --create-dirs -o ./datasets/dstress-16-plas.dat
python epnn.py

# linux
wget -c https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstate-16-plas.dat -P ./datasets/
wget -c https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstress-16-plas.dat -P ./datasets/
# windows
# curl https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstate-16-plas.dat --create-dirs -o ./datasets/dstate-16-plas.dat
# curl https://paddle-org.bj.bcebos.com/paddlescience/datasets/epnn/dstress-16-plas.dat --create-dirs -o ./datasets/dstress-16-plas.dat
python epnn.py mode=eval EVAL.pretrained_model_path=https://paddle-org.bj.bcebos.com/paddlescience/models/epnn/epnn_pretrained.pdparams
Pretrained Model Metrics
epnn_pretrained.pdparams error(total): 3.96903
error(error_elasto): 0.65328
error(error_plastic): 3.04176
error(error_stress): 0.27399

1. Background Introduction

Here we mainly reproduce the Physics-Informed Neural Network (PINN) surrogate model of the Elasto-Plastic Neural Network (EPNN). Incorporating these physics into the architecture of neural networks can train the network more effectively while using less data for training, and at the same time enhance inference capabilities for loading regimes outside the training data. The architecture of EPNN is model and material agnostic, meaning it can adapt to various types of elastoplastic materials, including geomaterials and metals; and experimental data can be used directly to train the network. To demonstrate the robustness of the proposed architecture, we apply its general framework to the elastoplastic behavior of sand. EPNN outperforms conventional neural network architectures in predicting unobserved strain-controlled loading paths for sands of different initial densities.

2. Problem Definition

In a neural network, information flows through connected neurons. The "strength" of each link in a neural network is determined by a variable weight:

\[ z_l^{\mathrm{i}}=W_{k l}^{\mathrm{i}-1, \mathrm{i}} a_k^{\mathrm{i}-1}+b^{\mathrm{i}-1}, \quad k=1: N^{\mathrm{i}-1} \quad \text { or } \quad \mathbf{z}^{\mathrm{i}}=\mathbf{a}^{\mathrm{i}-1} \mathbf{W}^{\mathrm{i}-1, \mathrm{i}}+b^{\mathrm{i}-1} \mathbf{I} \]

Where \(b\) is the bias term; \(N\) is the number of neurons in different layers; \(I\) refers to the unit vector where all elements are 1.

3. Problem Solving

Next, we will explain how to convert the problem into PaddleScience code step by step and solve the problem using deep learning methods. In order to quickly understand PaddleScience, only key steps such as model construction, equation construction, and computational domain construction are described below, while other details please refer to API Documentation.

3.1 Model Construction

In the EPNN problem, build the network, expressed in PaddleScience code as follows

371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
node_sizes_state_plastic.extend(hl_nodes_plastic)
node_sizes_stress.extend(hl_nodes_elasto)
node_sizes_state_elasto.extend([state_y_output_size - 3])
node_sizes_state_plastic.extend([state_y_output_size - 1])
node_sizes_stress.extend([1])

activation_elasto = "leaky_relu"
activation_plastic = "leaky_relu"
activations_elasto = [activation_elasto]
activations_plastic = [activation_plastic]
activations_elasto.extend([activation_elasto for ii in range(nhlayers)])
activations_plastic.extend([activation_plastic for ii in range(NHLAYERS_PLASTIC)])
activations_elasto.extend([activation_elasto])
activations_plastic.extend([activation_plastic])
drop_p = 0.0
n_state_elasto = ppsci.arch.Epnn(
    ("state_x",),
    ("out_state_elasto",),
    tuple(node_sizes_state_elasto),
    tuple(activations_elasto),
    drop_p,

EPNN parameters input_keys are input field names, output_keys are output field names, node_sizes are node size lists, activations are activation function string lists, and drop_p is node dropout probability.

3.2 Data Generation

This case involves reading data generation, as shown below

36
37
38
39
40
41
(
    input_dict_train,
    label_dict_train,
    input_dict_val,
    label_dict_val,
) = functions.get_data(cfg.DATASET_STATE, cfg.DATASET_STRESS, cfg.NTRAIN_SIZE)

306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
        n_train = math.floor(self.train_p * self.x.shape[0])
        n_cross_valid = math.floor(self.cross_valid_p * self.x.shape[0])
        n_test = math.floor(self.test_p * self.x.shape[0])
        self.x_train = self.x[shuffled_indices[0:n_train]]
        self.y_train = self.y[shuffled_indices[0:n_train]]
        self.x_valid = self.x[shuffled_indices[n_train : n_train + n_cross_valid]]
        self.y_valid = self.y[shuffled_indices[n_train : n_train + n_cross_valid]]
        self.x_test = self.x[
            shuffled_indices[n_train + n_cross_valid : n_train + n_cross_valid + n_test]
        ]
        self.y_test = self.y[
            shuffled_indices[n_train + n_cross_valid : n_train + n_cross_valid + n_test]
        ]


def get_data(dataset_state, dataset_stress, ntrain_size):

Here, Data is used to read files to construct the data class, then get_shuffled_data is used to shuffle the data, then the number of shuffled data itrain to be obtained is calculated, and finally get is used to obtain 10 groups of data with the quantity of itrain for each group.

3.3 Constraint Construction

Set training dataset and loss calculation function, return fields, code is as follows:

63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
output_keys = [
    "state_x",
    "state_y",
    "stress_x",
    "stress_y",
    "out_state_elasto",
    "out_state_plastic",
    "out_stress",
]
sup_constraint_pde = ppsci.constraint.SupervisedConstraint(
    {
        "dataset": {
            "name": "NamedArrayDataset",
            "input": input_dict_train,
            "label": label_dict_train,
        },
        "batch_size": 1,
        "num_workers": 0,
    },
    ppsci.loss.FunctionalLoss(functions.train_loss_func),
    {key: (lambda out, k=key: out[k]) for key in output_keys},
    name="sup_train",
)
constraint_pde = {sup_constraint_pde.name: sup_constraint_pde}

The first parameter of SupervisedConstraint is the reading configuration of the supervised constraint. The "dataset" field in the configuration represents the training dataset information used, and its various fields represent:

  1. name: Dataset type, here "NamedArrayDataset" means sequentially read dataset;
  2. input: Input dataset;
  3. label: Label dataset;

The second parameter is the loss function, here the custom function train_loss_func is used.

The third parameter is the equation expression, used to describe how to calculate the constraint target. The calculated value will be stored in the output list according to the specified name, so as to ensure that these values can be used when calculating loss.

The fourth parameter is the name of the constraint condition. We need to name each constraint condition for subsequent indexing.

After the constraint is constructed, encapsulate it into a dictionary with the name we just named as the key for subsequent access.

3.4 Validator Construction

Similar to constraints, this problem uses ppsci.validate.SupervisedValidator to build a validator. The parameter meanings are also similar to Constraint Construction. The only difference is the evaluation metric metric. The code is as follows:

 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
sup_validator_pde = ppsci.validate.SupervisedValidator(
    {
        "dataset": {
            "name": "NamedArrayDataset",
            "input": input_dict_val,
            "label": label_dict_val,
        },
        "batch_size": 1,
        "num_workers": 0,
    },
    ppsci.loss.FunctionalLoss(functions.eval_loss_func),
    {key: (lambda out, k=key: out[k]) for key in output_keys},
    metric={"metric": ppsci.metric.FunctionalMetric(functions.metric_expr)},
    name="sup_valid",
)
validator_pde = {sup_validator_pde.name: sup_validator_pde}

3.5 Hyperparameter Setting

Next we need to specify the number of training epochs. Here we use 10000 training epochs based on experimental experience. iters_per_epoch is 1.

40
41
TRAIN:
  epochs: 10000

3.6 Optimizer Construction

The training process will call the optimizer to update model parameters. Here, the more commonly used Adam optimizer is selected, and combined with the ExponentialDecay learning rate adjustment strategy commonly used in machine learning.

Since multiple models are used, multiple optimizers need to be set. For the EPNN network part, Adam optimizer needs to be set.

395
396
397
398
399
400
401
402
403
404
    ("out_state_plastic",),
    tuple(node_sizes_state_plastic),
    tuple(activations_plastic),
    drop_p,
)
n_stress = ppsci.arch.Epnn(
    ("state_x_f",),
    ("out_stress",),
    tuple(node_sizes_stress),
    tuple(activations_elasto),

Then for the added gkratio parameter, another optimizer needs to be set.

406
407
408
409
410
411
412
413
    )
    return (n_state_elasto, n_state_plastic, n_stress)


def get_optimizer_list(model_list, cfg):
    optimizer_list = []
    lr_list = [0.001, 0.001, 0.01]
    for i, model in enumerate(model_list):

Optimizers optimize in order, code summarized as:

395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
        ("out_state_plastic",),
        tuple(node_sizes_state_plastic),
        tuple(activations_plastic),
        drop_p,
    )
    n_stress = ppsci.arch.Epnn(
        ("state_x_f",),
        ("out_stress",),
        tuple(node_sizes_stress),
        tuple(activations_elasto),
        drop_p,
    )
    return (n_state_elasto, n_state_plastic, n_stress)


def get_optimizer_list(model_list, cfg):
    optimizer_list = []
    lr_list = [0.001, 0.001, 0.01]
    for i, model in enumerate(model_list):

3.7 Custom loss

Since this problem includes unsupervised learning and there is no label data in the data, the loss is calculated based on the returned data of the model, so a custom loss is required. The method is to first define relevant functions, and then pass the function name as a parameter to FunctionalLoss and FunctionalMetric.

Note that the input and output parameters of the custom loss function need to be consistent with other functions such as MSE in PaddleScience, that is, the input is the model output output_dict and other dictionary variables, and the loss function output is the loss value paddle.Tensor.

The relevant custom loss function is calculated using MAELoss, code is

114
115
116
117
118
119
120
121
122
123
124
125
126
        )
    }


def train_loss_func(output_dict, *args) -> paddle.Tensor:
    """For model calculation of loss in model.train().

    Args:
        output_dict (Dict[str, paddle.Tensor]): The output dict.

    Returns:
        paddle.Tensor: Loss value.
    """

3.8 Model Training and Evaluation

After completing the above settings, you only need to pass the instantiated objects to ppsci.solver.Solver in order.

106
107
108
109
110
111
112
113
114
115
116
117
118
solver = ppsci.solver.Solver(
    model_list_obj,
    constraint_pde,
    cfg.output_dir,
    optimizer_list,
    None,
    cfg.TRAIN.epochs,
    cfg.TRAIN.iters_per_epoch,
    save_freq=cfg.TRAIN.save_freq,
    eval_during_train=cfg.TRAIN.eval_during_train,
    validator=validator_pde,
    eval_with_no_grad=cfg.EVAL.eval_with_no_grad,
)

Set eval_during_train to True during model training, and evaluation will be performed after each training.

43
save_freq: 50

Finally start training:

121
solver.train()

4. Complete Code

epnn.py
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
# Copyright (c) 2023 PaddlePaddle Authors. All Rights Reserved.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""
Reference: https://github.com/meghbali/ANNElastoplasticity
"""

from os import path as osp

import functions
import hydra
from omegaconf import DictConfig

import ppsci
from ppsci.utils import logger


def train(cfg: DictConfig):
    # set random seed for reproducibility
    ppsci.utils.misc.set_random_seed(cfg.seed)

    # initialize logger
    logger.init_logger("ppsci", osp.join(cfg.output_dir, f"{cfg.mode}.log"), "info")

    (
        input_dict_train,
        label_dict_train,
        input_dict_val,
        label_dict_val,
    ) = functions.get_data(cfg.DATASET_STATE, cfg.DATASET_STRESS, cfg.NTRAIN_SIZE)
    model_list = functions.get_model_list(
        cfg.MODEL.ihlayers,
        cfg.MODEL.ineurons,
        input_dict_train["state_x"][0].shape[1],
        input_dict_train["state_y"][0].shape[1],
        input_dict_train["stress_x"][0].shape[1],
    )
    optimizer_list = functions.get_optimizer_list(model_list, cfg)
    model_state_elasto, model_state_plastic, model_stress = model_list
    model_list_obj = ppsci.arch.ModelList(model_list)

    def _transform_in_stress(_in):
        return functions.transform_in_stress(
            _in, model_state_elasto, "out_state_elasto"
        )

    model_state_elasto.register_input_transform(functions.transform_in)
    model_state_plastic.register_input_transform(functions.transform_in)
    model_stress.register_input_transform(_transform_in_stress)
    model_stress.register_output_transform(functions.transform_out)

    output_keys = [
        "state_x",
        "state_y",
        "stress_x",
        "stress_y",
        "out_state_elasto",
        "out_state_plastic",
        "out_stress",
    ]
    sup_constraint_pde = ppsci.constraint.SupervisedConstraint(
        {
            "dataset": {
                "name": "NamedArrayDataset",
                "input": input_dict_train,
                "label": label_dict_train,
            },
            "batch_size": 1,
            "num_workers": 0,
        },
        ppsci.loss.FunctionalLoss(functions.train_loss_func),
        {key: (lambda out, k=key: out[k]) for key in output_keys},
        name="sup_train",
    )
    constraint_pde = {sup_constraint_pde.name: sup_constraint_pde}

    sup_validator_pde = ppsci.validate.SupervisedValidator(
        {
            "dataset": {
                "name": "NamedArrayDataset",
                "input": input_dict_val,
                "label": label_dict_val,
            },
            "batch_size": 1,
            "num_workers": 0,
        },
        ppsci.loss.FunctionalLoss(functions.eval_loss_func),
        {key: (lambda out, k=key: out[k]) for key in output_keys},
        metric={"metric": ppsci.metric.FunctionalMetric(functions.metric_expr)},
        name="sup_valid",
    )
    validator_pde = {sup_validator_pde.name: sup_validator_pde}

    # initialize solver
    solver = ppsci.solver.Solver(
        model_list_obj,
        constraint_pde,
        cfg.output_dir,
        optimizer_list,
        None,
        cfg.TRAIN.epochs,
        cfg.TRAIN.iters_per_epoch,
        save_freq=cfg.TRAIN.save_freq,
        eval_during_train=cfg.TRAIN.eval_during_train,
        validator=validator_pde,
        eval_with_no_grad=cfg.EVAL.eval_with_no_grad,
    )

    # train model
    solver.train()
    functions.plotting(cfg.output_dir)


def evaluate(cfg: DictConfig):
    # set random seed for reproducibility
    ppsci.utils.misc.set_random_seed(cfg.seed)
    # initialize logger
    logger.init_logger("ppsci", osp.join(cfg.output_dir, f"{cfg.mode}.log"), "info")

    (
        input_dict_train,
        _,
        input_dict_val,
        label_dict_val,
    ) = functions.get_data(cfg.DATASET_STATE, cfg.DATASET_STRESS, cfg.NTRAIN_SIZE)
    model_list = functions.get_model_list(
        cfg.MODEL.ihlayers,
        cfg.MODEL.ineurons,
        input_dict_train["state_x"][0].shape[1],
        input_dict_train["state_y"][0].shape[1],
        input_dict_train["stress_x"][0].shape[1],
    )
    model_state_elasto, model_state_plastic, model_stress = model_list
    model_list_obj = ppsci.arch.ModelList(model_list)

    def _transform_in_stress(_in):
        return functions.transform_in_stress(
            _in, model_state_elasto, "out_state_elasto"
        )

    model_state_elasto.register_input_transform(functions.transform_in)
    model_state_plastic.register_input_transform(functions.transform_in)
    model_stress.register_input_transform(_transform_in_stress)
    model_stress.register_output_transform(functions.transform_out)

    output_keys = [
        "state_x",
        "state_y",
        "stress_x",
        "stress_y",
        "out_state_elasto",
        "out_state_plastic",
        "out_stress",
    ]
    sup_validator_pde = ppsci.validate.SupervisedValidator(
        {
            "dataset": {
                "name": "NamedArrayDataset",
                "input": input_dict_val,
                "label": label_dict_val,
            },
            "batch_size": 1,
            "num_workers": 0,
        },
        ppsci.loss.FunctionalLoss(functions.eval_loss_func),
        {key: (lambda out, k=key: out[k]) for key in output_keys},
        metric={"metric": ppsci.metric.FunctionalMetric(functions.metric_expr)},
        name="sup_valid",
    )
    validator_pde = {sup_validator_pde.name: sup_validator_pde}
    functions.OUTPUT_DIR = cfg.output_dir

    # initialize solver
    solver = ppsci.solver.Solver(
        model_list_obj,
        output_dir=cfg.output_dir,
        validator=validator_pde,
        pretrained_model_path=cfg.EVAL.pretrained_model_path,
        eval_with_no_grad=cfg.EVAL.eval_with_no_grad,
    )
    # evaluate
    solver.eval()


@hydra.main(version_base=None, config_path="./conf", config_name="epnn.yaml")
def main(cfg: DictConfig):
    if cfg.mode == "train":
        train(cfg)
    elif cfg.mode == "eval":
        evaluate(cfg)
    else:
        raise ValueError(f"cfg.mode should in ['train', 'eval'], but got '{cfg.mode}'")


if __name__ == "__main__":
    main()

5. Result Display

The EPNN case was experimented with the parameter configuration of epoch=10000, and the result returned Loss was 0.00471.

The figures below are the Loss, Training error, and Cross validation error graphs for different epochs:

loss_trend

Training loss graph

6. References