Skip to content

2D-Darcy

AI Studio Quick Experience

python darcy2d.py

python darcy2d.py mode=eval EVAL.pretrained_model_path=https://paddle-org.bj.bcebos.com/paddlescience/models/darcy2d/darcy2d_pretrained.pdparams

python darcy2d.py mode=export

python darcy2d.py mode=infer
Pretrained Model Metrics
darcy2d_pretrained.pdparams loss(Residual): 0.36500
MSE.poisson(Residual): 0.00006

1. Background Introduction

Darcy flow describes fluid movement through porous media based on Darcy's Law, finding extensive applications in groundwater modeling, hydrology, hydrogeology, and petroleum engineering.

In petroleum engineering, Darcy flow is essential for predicting and simulating oil movement within porous rock formations. By modeling flow through the voids between particles, engineers can optimize extraction and production processes.

In hydrogeology and agriculture, Darcy flow models are used to study groundwater dynamics. This includes predicting the effects of irrigation on soil moisture to optimize crop planning, as well as forecasting and mitigating groundwater pollution in urban planning.

The 2D-Darcy problem focuses on linear seepage in a two-dimensional plane. When fluid velocity in porous media is low, flow adheres to Darcy's law, exhibiting a linear relationship between seepage velocity and the pressure gradient.

2. Problem Definition

Assume that in the Darcy flow model, the flow velocity \(\mathbf{u}\) and pressure \(p\) at each position \((x,y)\) satisfy the following relationship:

\[ \begin{cases} \begin{aligned} \mathbf{u}+\nabla p =& 0,(x,y) \in \Omega \\ \nabla \cdot \mathbf{u} =& f,(x,y) \in \Omega \\ p(x,y) =& \sin(2 \pi x )\cos(2 \pi y), (x,y) \in \partial \Omega \end{aligned} \end{cases} \]

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 2D-Darcy problem, for each coordinate \((x, y)\), the unknown pressure \(p\) must be determined. We use a Multilayer Perceptron (MLP) to approximate the mapping function \(f: \mathbb{R}^2 \to \mathbb{R}^1\) from \((x, y)\) to \(p\), such that:

\[ p = f(x, y) \]

In the above formula, \(f\) is the MLP model itself, expressed in PaddleScience code as follows:

33
34
# set model
model = ppsci.arch.MLP(**cfg.MODEL)

To ensure accurate and efficient variable access during computation, we define the model's input keys as ("x", "y") and the output key as "p", maintaining consistency with the code.

We instantiate the MLP model with 5 hidden layers and 20 neurons per layer.

3.2 Equation Construction

Since this problem involves the 2D Poisson equation, we can directly use the built-in Poisson class in PaddleScience, setting the dim parameter to 2.

36
37
# set equation
equation = {"Poisson": ppsci.equation.Poisson(2)}

3.3 Computational Domain Construction

The 2D-Darcy problem is defined on a rectangular domain with diagonal corners at (0.0, 0.0) and (1.0, 1.0). We can directly use the built-in Rectangle geometry in PaddleScience.

39
40
# set geometry
geom = {"rect": ppsci.geometry.Rectangle((0.0, 0.0), (1.0, 1.0))}

3.4 Constraint Construction

We use two constraints to guide model training: the Darcy equation constraint on internal points and boundary constraints on the domain edges.

Before defining these constraints, we specify the sampling configuration, including the number of points for each constraint.

42
43
44
45
46
# set dataloader config
train_dataloader_cfg = {
    "dataset": "IterableNamedArrayDataset",
    "iters_per_epoch": cfg.TRAIN.iters_per_epoch,
}

3.4.1 Interior Point Constraint

Taking InteriorConstraint acting on internal points as an example, the code is as follows:

48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
# set constraint
def poisson_ref_compute_func(_in):
    return (
        -8.0
        * (np.pi**2)
        * np.sin(2.0 * np.pi * _in["x"])
        * np.cos(2.0 * np.pi * _in["y"])
    )

pde_constraint = ppsci.constraint.InteriorConstraint(
    equation["Poisson"].equations,
    {"poisson": poisson_ref_compute_func},
    geom["rect"],
    {**train_dataloader_cfg, "batch_size": cfg.NPOINT_PDE},
    ppsci.loss.MSELoss("sum"),
    evenly=True,
    name="EQ",
)
  • Equation Expression: Specifies how to calculate the constraint target. We use equation["Poisson"].equations from Section 3.2.
  • Target Values: The target values for the constraint variables. We set these to the results generated by poisson_ref_compute_func to match the standard solution.
  • Computational Domain: The domain where the constraint applies. We use geom["rect"] from Section 3.3.
  • Sampling Configuration: We use full data points (IterableNamedArrayDataset, iters_per_epoch=1) with a batch_size of 9801 (a 99x99 grid).
  • Loss Function: We use the MSE function with reduction="sum" to sum the loss across all points.
  • Equidistant Sampling: We enable equidistant sampling to distribute training points evenly across the domain, facilitating convergence.
  • Name: A unique name for the constraint, set here as "EQ".

3.4.2 Boundary Constraint

We also construct constraints for the four boundaries of the rectangle. Since these apply to the boundary, we use the BoundaryConstraint class:

67
68
69
70
71
72
73
74
75
76
77
bc = ppsci.constraint.BoundaryConstraint(
    {"p": lambda out: out["p"]},
    {
        "p": lambda _in: np.sin(2.0 * np.pi * _in["x"])
        * np.cos(2.0 * np.pi * _in["y"])
    },
    geom["rect"],
    {**train_dataloader_cfg, "batch_size": cfg.NPOINT_BC},
    ppsci.loss.MSELoss("sum"),
    name="BC",
)

The first parameter specifies the constraint object, which is the model output out["p"].

The second parameter defines the target value, calculated directly via the analytical solution:

lambda _in: np.sin(2.0 * np.pi * _in["x"]) * np.cos(2.0 * np.pi * _in["y"])

Other parameters for BoundaryConstraint are similar to those in InteriorConstraint.

Finally, we encapsulate the constraints into a dictionary for easy access.

78
79
80
81
82
# wrap constraints together
constraint = {
    pde_constraint.name: pde_constraint,
    bc.name: bc,
}

3.5 Hyperparameter Setting

Next, we need to specify the number of training epochs and learning rate. Here, based on experimental experience, we use 10,000 training epochs.

39
40
41
42
43
44
45
46
47
  activation: "stan"

# training settings
TRAIN:
  epochs: 10000
  iters_per_epoch: 1
  lr_scheduler:
    epochs: ${TRAIN.epochs}
    iters_per_epoch: ${TRAIN.iters_per_epoch}

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 the OneCycle learning rate adjustment strategy commonly used in machine learning is used together.

84
85
86
# set optimizer
lr_scheduler = ppsci.optimizer.lr_scheduler.OneCycleLR(**cfg.TRAIN.lr_scheduler)()
optimizer = ppsci.optimizer.Adam(lr_scheduler)(model)

3.7 Validator Construction

Usually during the training process, the training status of the current model is evaluated using the validation set (test set) at a certain epoch interval, so ppsci.validate.GeometryValidator is used to construct the validator.

 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
# set validator
residual_validator = ppsci.validate.GeometryValidator(
    equation["Poisson"].equations,
    {"poisson": poisson_ref_compute_func},
    geom["rect"],
    {
        "dataset": "NamedArrayDataset",
        "total_size": cfg.NPOINT_PDE,
        "batch_size": cfg.EVAL.batch_size.residual_validator,
        "sampler": {"name": "BatchSampler"},
    },
    ppsci.loss.MSELoss("sum"),
    evenly=True,
    metric={"MSE": ppsci.metric.MSE()},
    name="Residual",
)
validator = {residual_validator.name: residual_validator}

3.8 Visualizer Construction

During model evaluation, if the evaluation result is data that can be visualized, we can choose a suitable visualizer to visualize the output result.

The output data in this article is a two-dimensional point set in an area, so we only need to save the evaluated output data as a vtu format file, and finally open it with visualization software to view. The code is as follows:

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
# set visualizer(optional)
# manually collate input data for visualization,
vis_points = geom["rect"].sample_interior(
    cfg.NPOINT_PDE + cfg.NPOINT_BC, evenly=True
)
visualizer = {
    "visualize_p_ux_uy": ppsci.visualize.VisualizerVtu(
        vis_points,
        {
            "p": lambda d: d["p"],
            "p_ref": lambda d: paddle.sin(2 * np.pi * d["x"])
            * paddle.cos(2 * np.pi * d["y"]),
            "p_diff": lambda d: paddle.sin(2 * np.pi * d["x"])
            * paddle.cos(2 * np.pi * d["y"])
            - d["p"],
            "ux": lambda d: jacobian(d["p"], d["x"]),
            "ux_ref": lambda d: 2
            * np.pi
            * paddle.cos(2 * np.pi * d["x"])
            * paddle.cos(2 * np.pi * d["y"]),
            "ux_diff": lambda d: jacobian(d["p"], d["x"])
            - 2
            * np.pi
            * paddle.cos(2 * np.pi * d["x"])
            * paddle.cos(2 * np.pi * d["y"]),
            "uy": lambda d: jacobian(d["p"], d["y"]),
            "uy_ref": lambda d: -2
            * np.pi
            * paddle.sin(2 * np.pi * d["x"])
            * paddle.sin(2 * np.pi * d["y"]),
            "uy_diff": lambda d: jacobian(d["p"], d["y"])
            - (
                -2
                * np.pi
                * paddle.sin(2 * np.pi * d["x"])
                * paddle.sin(2 * np.pi * d["y"])
            ),
        },
        prefix="result_p_ux_uy",
    )
}

3.9 Model Training, Evaluation and Visualization

3.9.1 Training with Adam

After completing the above settings, you only need to pass the instantiated objects to ppsci.solver.Solver in order, and then start training, evaluation, and visualization.

148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
# initialize solver
solver = ppsci.solver.Solver(
    model,
    constraint,
    cfg.output_dir,
    optimizer,
    lr_scheduler,
    cfg.TRAIN.epochs,
    cfg.TRAIN.iters_per_epoch,
    eval_during_train=cfg.TRAIN.eval_during_train,
    eval_freq=cfg.TRAIN.eval_freq,
    equation=equation,
    geom=geom,
    validator=validator,
    visualizer=visualizer,
)
# train model
solver.train()
# evaluate after finished training
solver.eval()
# visualize prediction after finished training
solver.visualize()

3.9.2 Finetuning with L-BFGS [Optional]

After training with the Adam optimizer, we can replace the optimizer with the second-order optimizer L-BFGS to continue training for a small number of rounds (here we use 10% of the Adam optimization rounds) to further improve model accuracy.

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
OUTPUT_DIR = cfg.TRAIN.lbfgs.output_dir
logger.init_logger("ppsci", osp.join(OUTPUT_DIR, f"{cfg.mode}.log"), "info")
EPOCHS = cfg.TRAIN.epochs // 10
optimizer_lbfgs = ppsci.optimizer.LBFGS(
    cfg.TRAIN.lbfgs.learning_rate, cfg.TRAIN.lbfgs.max_iter
)(model)
solver = ppsci.solver.Solver(
    model,
    constraint,
    OUTPUT_DIR,
    optimizer_lbfgs,
    None,
    EPOCHS,
    cfg.TRAIN.lbfgs.iters_per_epoch,
    eval_during_train=cfg.TRAIN.lbfgs.eval_during_train,
    eval_freq=cfg.TRAIN.lbfgs.eval_freq,
    equation=equation,
    geom=geom,
    validator=validator,
    visualizer=visualizer,
)
# train model
solver.train()
# evaluate after finished training
solver.eval()
# visualize prediction after finished training
solver.visualize()
Tip

Using L-BFGS to fine-tune for a small number of rounds after training with conventional optimizers can further effectively improve model accuracy in most scenarios.

4. Complete Code

darcy2d.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
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
# 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.

from os import path as osp

import hydra
import numpy as np
import paddle
from omegaconf import DictConfig

import ppsci
from ppsci.autodiff import jacobian
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")

    # set model
    model = ppsci.arch.MLP(**cfg.MODEL)

    # set equation
    equation = {"Poisson": ppsci.equation.Poisson(2)}

    # set geometry
    geom = {"rect": ppsci.geometry.Rectangle((0.0, 0.0), (1.0, 1.0))}

    # set dataloader config
    train_dataloader_cfg = {
        "dataset": "IterableNamedArrayDataset",
        "iters_per_epoch": cfg.TRAIN.iters_per_epoch,
    }

    # set constraint
    def poisson_ref_compute_func(_in):
        return (
            -8.0
            * (np.pi**2)
            * np.sin(2.0 * np.pi * _in["x"])
            * np.cos(2.0 * np.pi * _in["y"])
        )

    pde_constraint = ppsci.constraint.InteriorConstraint(
        equation["Poisson"].equations,
        {"poisson": poisson_ref_compute_func},
        geom["rect"],
        {**train_dataloader_cfg, "batch_size": cfg.NPOINT_PDE},
        ppsci.loss.MSELoss("sum"),
        evenly=True,
        name="EQ",
    )

    bc = ppsci.constraint.BoundaryConstraint(
        {"p": lambda out: out["p"]},
        {
            "p": lambda _in: np.sin(2.0 * np.pi * _in["x"])
            * np.cos(2.0 * np.pi * _in["y"])
        },
        geom["rect"],
        {**train_dataloader_cfg, "batch_size": cfg.NPOINT_BC},
        ppsci.loss.MSELoss("sum"),
        name="BC",
    )
    # wrap constraints together
    constraint = {
        pde_constraint.name: pde_constraint,
        bc.name: bc,
    }

    # set optimizer
    lr_scheduler = ppsci.optimizer.lr_scheduler.OneCycleLR(**cfg.TRAIN.lr_scheduler)()
    optimizer = ppsci.optimizer.Adam(lr_scheduler)(model)

    # set validator
    residual_validator = ppsci.validate.GeometryValidator(
        equation["Poisson"].equations,
        {"poisson": poisson_ref_compute_func},
        geom["rect"],
        {
            "dataset": "NamedArrayDataset",
            "total_size": cfg.NPOINT_PDE,
            "batch_size": cfg.EVAL.batch_size.residual_validator,
            "sampler": {"name": "BatchSampler"},
        },
        ppsci.loss.MSELoss("sum"),
        evenly=True,
        metric={"MSE": ppsci.metric.MSE()},
        name="Residual",
    )
    validator = {residual_validator.name: residual_validator}

    # set visualizer(optional)
    # manually collate input data for visualization,
    vis_points = geom["rect"].sample_interior(
        cfg.NPOINT_PDE + cfg.NPOINT_BC, evenly=True
    )
    visualizer = {
        "visualize_p_ux_uy": ppsci.visualize.VisualizerVtu(
            vis_points,
            {
                "p": lambda d: d["p"],
                "p_ref": lambda d: paddle.sin(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "p_diff": lambda d: paddle.sin(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"])
                - d["p"],
                "ux": lambda d: jacobian(d["p"], d["x"]),
                "ux_ref": lambda d: 2
                * np.pi
                * paddle.cos(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "ux_diff": lambda d: jacobian(d["p"], d["x"])
                - 2
                * np.pi
                * paddle.cos(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "uy": lambda d: jacobian(d["p"], d["y"]),
                "uy_ref": lambda d: -2
                * np.pi
                * paddle.sin(2 * np.pi * d["x"])
                * paddle.sin(2 * np.pi * d["y"]),
                "uy_diff": lambda d: jacobian(d["p"], d["y"])
                - (
                    -2
                    * np.pi
                    * paddle.sin(2 * np.pi * d["x"])
                    * paddle.sin(2 * np.pi * d["y"])
                ),
            },
            prefix="result_p_ux_uy",
        )
    }

    # initialize solver
    solver = ppsci.solver.Solver(
        model,
        constraint,
        cfg.output_dir,
        optimizer,
        lr_scheduler,
        cfg.TRAIN.epochs,
        cfg.TRAIN.iters_per_epoch,
        eval_during_train=cfg.TRAIN.eval_during_train,
        eval_freq=cfg.TRAIN.eval_freq,
        equation=equation,
        geom=geom,
        validator=validator,
        visualizer=visualizer,
    )
    # train model
    solver.train()
    # evaluate after finished training
    solver.eval()
    # visualize prediction after finished training
    solver.visualize()

    # fine-tuning pretrained model with L-BFGS
    OUTPUT_DIR = cfg.TRAIN.lbfgs.output_dir
    logger.init_logger("ppsci", osp.join(OUTPUT_DIR, f"{cfg.mode}.log"), "info")
    EPOCHS = cfg.TRAIN.epochs // 10
    optimizer_lbfgs = ppsci.optimizer.LBFGS(
        cfg.TRAIN.lbfgs.learning_rate, cfg.TRAIN.lbfgs.max_iter
    )(model)
    solver = ppsci.solver.Solver(
        model,
        constraint,
        OUTPUT_DIR,
        optimizer_lbfgs,
        None,
        EPOCHS,
        cfg.TRAIN.lbfgs.iters_per_epoch,
        eval_during_train=cfg.TRAIN.lbfgs.eval_during_train,
        eval_freq=cfg.TRAIN.lbfgs.eval_freq,
        equation=equation,
        geom=geom,
        validator=validator,
        visualizer=visualizer,
    )
    # train model
    solver.train()
    # evaluate after finished training
    solver.eval()
    # visualize prediction after finished training
    solver.visualize()


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")

    # set model
    model = ppsci.arch.MLP(**cfg.MODEL)

    # set equation
    equation = {"Poisson": ppsci.equation.Poisson(2)}

    # set geometry
    geom = {"rect": ppsci.geometry.Rectangle((0.0, 0.0), (1.0, 1.0))}

    # set constraint
    def poisson_ref_compute_func(_in):
        return (
            -8.0
            * (np.pi**2)
            * np.sin(2.0 * np.pi * _in["x"])
            * np.cos(2.0 * np.pi * _in["y"])
        )

    # set validator
    residual_validator = ppsci.validate.GeometryValidator(
        equation["Poisson"].equations,
        {"poisson": poisson_ref_compute_func},
        geom["rect"],
        {
            "dataset": "NamedArrayDataset",
            "total_size": cfg.NPOINT_PDE,
            "batch_size": cfg.EVAL.batch_size.residual_validator,
            "sampler": {"name": "BatchSampler"},
        },
        ppsci.loss.MSELoss("sum"),
        evenly=True,
        metric={"MSE": ppsci.metric.MSE()},
        name="Residual",
    )
    validator = {residual_validator.name: residual_validator}

    # set visualizer
    # manually collate input data for visualization,
    vis_points = geom["rect"].sample_interior(
        cfg.NPOINT_PDE + cfg.NPOINT_BC, evenly=True
    )
    visualizer = {
        "visualize_p_ux_uy": ppsci.visualize.VisualizerVtu(
            vis_points,
            {
                "p": lambda d: d["p"],
                "p_ref": lambda d: paddle.sin(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "p_diff": lambda d: paddle.sin(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"])
                - d["p"],
                "ux": lambda d: jacobian(d["p"], d["x"]),
                "ux_ref": lambda d: 2
                * np.pi
                * paddle.cos(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "ux_diff": lambda d: jacobian(d["p"], d["x"])
                - 2
                * np.pi
                * paddle.cos(2 * np.pi * d["x"])
                * paddle.cos(2 * np.pi * d["y"]),
                "uy": lambda d: jacobian(d["p"], d["y"]),
                "uy_ref": lambda d: -2
                * np.pi
                * paddle.sin(2 * np.pi * d["x"])
                * paddle.sin(2 * np.pi * d["y"]),
                "uy_diff": lambda d: jacobian(d["p"], d["y"])
                - (
                    -2
                    * np.pi
                    * paddle.sin(2 * np.pi * d["x"])
                    * paddle.sin(2 * np.pi * d["y"])
                ),
            },
            prefix="result_p_ux_uy",
        )
    }

    solver = ppsci.solver.Solver(
        model,
        output_dir=cfg.output_dir,
        equation=equation,
        geom=geom,
        validator=validator,
        visualizer=visualizer,
        pretrained_model_path=cfg.EVAL.pretrained_model_path,
    )
    solver.eval()
    # visualize prediction
    solver.visualize()


def export(cfg: DictConfig):
    # set model
    model = ppsci.arch.MLP(**cfg.MODEL)

    # initialize solver
    solver = ppsci.solver.Solver(
        model,
        pretrained_model_path=cfg.INFER.pretrained_model_path,
    )
    # export model
    from paddle.static import InputSpec

    input_spec = [
        {key: InputSpec([None, 1], "float32", name=key) for key in model.input_keys},
    ]

    solver.export(input_spec, cfg.INFER.export_path)


def inference(cfg: DictConfig):
    from deploy.python_infer import pinn_predictor

    predictor = pinn_predictor.PINNPredictor(cfg)

    # set geometry
    geom = {"rect": ppsci.geometry.Rectangle((0.0, 0.0), (1.0, 1.0))}
    # manually collate input data for visualization,
    input_dict = geom["rect"].sample_interior(
        cfg.NPOINT_PDE + cfg.NPOINT_BC, evenly=True
    )
    output_dict = predictor.predict(
        {key: input_dict[key] for key in cfg.MODEL.input_keys}, cfg.INFER.batch_size
    )
    # mapping data to cfg.INFER.output_keys
    output_dict = {
        store_key: output_dict[infer_key]
        for store_key, infer_key in zip(cfg.MODEL.output_keys, output_dict.keys())
    }
    ppsci.visualize.save_vtu_from_dict(
        "./visual/darcy2d.vtu",
        {**input_dict, **output_dict},
        input_dict.keys(),
        cfg.MODEL.output_keys,
    )


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


if __name__ == "__main__":
    main()

5. Result Display

The figures below display the prediction results, reference values, and errors for pressure \(p(x,y)\), x-velocity \(u(x,y)\), and y-velocity \(v(x,y)\) across the computational domain.

darcy 2d

Left: Predicted pressure p, Middle: Reference pressure p, Right: Pressure difference

darcy 2d

Left: Predicted x-direction velocity p, Middle: Reference x-direction velocity p, Right: x-direction velocity difference

darcy 2d

Left: Predicted y-direction velocity p, Middle: Reference y-direction velocity p, Right: y-direction velocity difference