&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&ensp;
[Home Page](../../Start_Here.ipynb)

[Previous Notebook](../introduction/Introductory_Notebook.ipynb)
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&ensp;
[1](../introduction/Introductory_Notebook.ipynb)
[2]
[3](../spring_mass/Spring_Mass_Problem_Notebook.ipynb)
[4](../chip_2d/Challenge_CFD_Problem_Notebook.ipynb)
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
[Next Notebook](../spring_mass/Spring_Mass_Problem_Notebook.ipynb)


# Steady State 1D Diffusion in a Composite Bar using PINNs

This notebook give you a headstart in solving your own Partial Differential Equations (PDEs) using neural networks. Let's quickly recap the theory of PINNs before proceeding. We will embed the physics of the problem in the the neural networks and in such setting, we can use them to approximate the solution to a given differential equation and boundary condition without any training data from other solvers. More specifically, the neural network will be trained to minimize a loss function which is formed using the differential equation and the boundary conditions. If the network is able to minimize this loss then it will in effect solve the given differential equation. More information about the Physics Informed Neural Networks (PINNs) can be found in the [paper](https://www.sciencedirect.com/science/article/pii/S0021999118307125?casa_token=1CnSbVeDwJ0AAAAA:-8a6ZMjO7RgYjFBAxIoVU2sWAdsqFzLRTNHA1BkRryNPXx5Vjc8hCDCkS99gcZR0B1rpa33a30yG) published by Raissi et al. 

In this notebook we will solve the steady 1 dimensional heat transfer in a composite bar. We will use NVIDIA's SimNet library to create the problem setup. You can refer to the *SimNet User Guide* for more examples on solving different types of PDEs using the SimNet library. Also, for more information about the SimNet APIs you can refer the *SimNet Source Code Documentation*. 

### Learning Outcomes
1. How to use SimNet to simulate physics problems using PINNs
    1. How to write your own PDEs and formulate the different losses
    2. How to use the Constructive Solid Geometry (CSG) module
2. How to use SimNet to solve a parameterized PDEs
    

## Problem Description

Our aim is to obtain the temperature distribution inside the bar that is made up of two materials with different thermal conductivity. The geometry and the problem specification of the problem can be seen below



The composite bar extends from $x=0$ to $x=2$. The bar has material of conductivity $D_1=10$ from $x=0$ to $x=1$ and $D_2=0.1$ from $x=1$ to $x=2$. Both the ends of the bar, $x=0$ and $x=2$ are maintained at a constant temperatures of $0$ and $100$ respectively. For simplicity of modeling, we will treat the composite bar as two separate bars, bar 1 and bar 2, whose ends are joined together. We will treat the temperatures in the bar 1 as $U_1$ and the temperature in bar 2 as $U_2$. 

The equations and boundary conditions governing the problem can be mathematically expressed as

One dimensional diffusion of temperature in bar 1 and 2:

$$
\begin{align}
\frac{d}{dx}\left( D_1\frac{dU_1}{dx} \right) = 0, && \text{when } 0<x<1 \\
\frac{d}{dx}\left( D_2\frac{dU_2}{dx} \right) = 0, && \text{when } 1<x<2 \\
\end{align}
$$

Flux and temperature continuity at interface $(x=1)$
$$
\begin{align}
D_1\frac{dU_1}{dx} = D_1\frac{dU_1}{dx}, && \text{when } x=1 \\
U_1 = U_2, && \text{when } x=1 \\
\end{align}
$$


## Case Setup

Now that we have our problem defined, let's take a look at the code required to solve it using SimNet's PINN library. SimNet has a variety of helper functions that will help us to set up the problem with ease. It has APIs to model geometry in a parameterized fashion using the Constructive Solid Geometry (CSG) module, write-up the required equations in a user-friendly symbolic format and comes with several advanced neural network architectures to choose for more complicated problems. 

Now let's start the problem by importing the required libraries and packages


```python
# import SimNet library
from sympy import Symbol, Function, Number, sin, Eq, Abs, exp
import numpy as np
import tensorflow as tf

from simnet.solver import Solver
from simnet.dataset import TrainDomain, ValidationDomain, MonitorDomain
from simnet.data import Validation, Monitor
from simnet.sympy_utils.geometry_1d import Line1D
from simnet.controller import SimNetController
from simnet.node import Node
from simnet.pdes import PDES
from simnet.variables import Variables
```

The `Solver` class trains and evaluates the SimNet's neural network solver. The class `TrainDomain` is used to define the training data for the problem, while the other classes like `ValidationDomain`, `MonitorDomain`, etc. are used to create other data evaluations during the training. 

The modules like `PDES` and `sympy_utils` contain predefined differential equations and geometries respectively that one can use to define the problem. We will describe each of them in detail as we move forward in the code. For more detailed information on all the different modules present in SimNet, we recommended you to refer the *SimNet Source Code Documentation*.



## Creating the geometry

In this problem, we will create the 1-dimensional geometry using the `Line1D` from the geometry module. The module also contains several 2d and 3d shapes like rectangle, circle, triangle, cuboids, sphere, torus, cones, tetrahedrons, etc. We will define the one dimensional line object using the two end-points. For composite bar, we will create two separate bars as defined in the problem statement


```python
# params for domain
L1 = Line1D(0,1)
L2 = Line1D(1,2)
```

Next we will define the properties for the problem which will later be used while making the equations, boundary conditions, etc. Also, for this problem, we can find the temperature at the interface analytically and we will use that as to form validation domain to compare our neural network results 


```python
# defining the parameters for boundary conditions and equations of the problem
D1 = 1e1
D2 = 1e-1
Ta = 0
Tc = 100

# temperature at the interface from analytical solution
Tb = (Tc + (D1/D2)*Ta)/(1 + (D1/D2))
```

## Defining the differential equations for the problem

The `PDES` class allows us to write the equations symbolically in Sympy. This allows users to quickly write their equations in the most natural way possible. The Sympy equations are converted to TensorFlow expressions in the back-end and can also be printed to ensure correct implementation.

SimNet also comes with several common PDEs predefined for the user to choose from. Some of the PDEs that are already available in the PDEs module are: Navier Stokes, Linear Elasticity, Advection Diffusion, Wave Equations, etc.

Let's create the PDE to define the diffusion equation. We will define the equation in its most generic, transient 3-dimensional, form and then have an argument `dim` that can reduce it to lower dimensional forms. 
$$\frac{\partial T}{\partial t}= \nabla\cdot \left( D \nabla T \right) + Q$$

Let's start defining the equation by inhereting from the `PDES` class. We will create the initialization method for this class that defines the equation(s) of interest. We will be defining the diffusion equation using the source(`Q`), diffusivity(`D`), symbol for diffusion(`T`). If `D` or `Q` is given as a string we will convert it to functional form. This will allow us to solve problems with spatially/temporally varying properties. 


```python
class Diffusion(PDES):
  name = 'Diffusion'
 
  def __init__(self, T='T', D='D', Q=0, dim=3, time=True):
    # set params
    self.T = T
    self.dim = dim
    self.time = time

    # coordinates
    x, y, z = Symbol('x'), Symbol('y'), Symbol('z')

    # time
    t = Symbol('t')

    # make input variables 
    input_variables = {'x':x,'y':y,'z':z,'t':t}
    if self.dim == 1:
      input_variables.pop('y')
      input_variables.pop('z')
    elif self.dim == 2:
      input_variables.pop('z')
    if not self.time: 
      input_variables.pop('t')

    # Temperature
    assert type(T) == str, "T needs to be string"
    T = Function(T)(*input_variables)

    # Diffusivity
    if type(D) is str:
      D = Function(D)(*input_variables)
    elif type(D) in [float, int]:
      D = Number(D)

    # Source
    if type(Q) is str:
      Q = Function(Q)(*input_variables)
    elif type(Q) in [float, int]:
      Q = Number(Q)

    # set equations
    self.equations = Variables()
    self.equations['diffusion_'+self.T] = (T.diff(t)
                                          - (D*T.diff(x)).diff(x)
                                          - (D*T.diff(y)).diff(y)
                                          - (D*T.diff(z)).diff(z)
                                          - Q)
```

First we defined the input variables $x, y, z$ and $t$ with Sympy symbols. Then we defined the functions for $T$, $D$ and $Q$ that are dependent on the input variables $(x, y, z, t)$. Using these we can write out our simple equation $T_t = \nabla \cdot (D \nabla T) + Q$. We store this equation in the class by adding it to the dictionary of `equations`.

Note that we moved all the terms of the PDE either to LHS or RHS. This way, while using the equations in the `TrainDomain`, we
can assign a custom source function to the `’diffusion_T’` key instead of 0 to add more source terms to our PDE.

Great! We just wrote our own PDE in SimNet! Once you have understood the process to code a simple PDE, you can easily extend the procedure for different PDEs. You can also bundle multiple PDEs together in a same file by adding new keys to the equations dictionary. Below we show the code for the interface boundary condition where we need to maintain the field (dirichlet) and flux (neumann) continuity. *(More examples of coding your own PDE can be found in the SimNet User Guide Chapter 4)*. 

**Note :** The field continuity condition is needed because we are solving for two different temperatures in the two bars. 


```python
class DiffusionInterface(PDES):
  name = 'DiffusionInterface'

  def __init__(self, T_1, T_2, D_1, D_2, dim=3, time=True):
    # set params
    self.T_1 = T_1
    self.T_2 = T_2
    self.D_1 = D_1
    self.D_2 = D_2
    self.dim = dim
    self.time = time
 
    # coordinates
    x, y, z = Symbol('x'), Symbol('y'), Symbol('z')
    normal_x, normal_y, normal_z = Symbol('normal_x'), Symbol('normal_y'), Symbol('normal_z')

    # time
    t = Symbol('t')

    # make input variables 
    input_variables = {'x':x,'y':y,'z':z,'t':t}
    if self.dim == 1:
      input_variables.pop('y')
      input_variables.pop('z')
    elif self.dim == 2:
      input_variables.pop('z')
    if not self.time: 
      input_variables.pop('t')

    # variables to match the boundary conditions (example Temperature)
    T_1 = Function(T_1)(*input_variables)
    T_2 = Function(T_2)(*input_variables)

    # set equations
    self.equations = Variables()
    self.equations['diffusion_interface_dirichlet_'+self.T_1+'_'+self.T_2] = T_1 - T_2
    flux_1 = self.D_1 * (normal_x * T_1.diff(x) + normal_y * T_1.diff(y) + normal_z * T_1.diff(z))
    flux_2 = self.D_2 * (normal_x * T_2.diff(x) + normal_y * T_2.diff(y) + normal_z * T_2.diff(z))
    self.equations['diffusion_interface_neumann_'+self.T_1+'_'+self.T_2] = flux_1 - flux_2
```

## Creating Train Domain: Assigning the boundary conditions and equations to the geometry

As described earlier, we need to define a training domain for training our neural network. A loss function is then constructed which is a combination of contributions from the boundary conditions and equations that a neural network must satisfy at the end of the training. These training points (BCs and equations) are defined in a class that inherits from the `TrainDomain` parent class. The boundary conditions are implemented as soft constraints. These BCs along with the equations to be solved are used to formulate a composite loss that is minimized by the network during training. 

$$L = L_{BC} + L_{Residual}$$

**Boundary conditions:** For generating a boundary condition, we need to sample the points on the required boundary/surface of the geometry and then assign them the desired values. We will use the method `boundary_bc` to sample the points on the boundary of the geometry we already created. `boundary_bc` will sample the entire boundary of the geometry, in this case, both the endpoints of the 1d line. A particular boundary of the geometry can be sub-sampled by using a particular criterion for the boundary_bc using the criteria parameter. For example, to sample the left end of `L1`, criteria is set to `Eq(x, 0)`. 

The desired values for the boundary condition are listed as a dictionary in `outvar_sympy` parameter. In SimNet we define these variables as keys of this dictionary which are converted to appropriate nodes in the computational graph. For this problem,  we have `'u_1':0` at $x=0$ and `'u_2':100` at $x=2$. At $x=1$, we have the interface condition `'diffusion_interface_dirichlet_u_1_u_2':0` and `'diffusion_interface_neumann_u_1_u_2':0` that we defined earlier (i.e. $U_1=U_2$ and $D_1\frac{dU_1}{dx}=D_2\frac{dU_2}{dx}$). 

The number of points to sample on each boundary are specified using the `batch_size_per_area` parameter. The
actual number of points sampled is then equal to the length/area of the geometry being sampled (boundary or interior)
times the batch_size_per_area. 

In this case, since we only have 1 point on the boundary, we specify the `batch_size_per_area` as 1 for all the boundaries.

**Equations to solve:** The Diffusion PDE we defined is enforced on all the points in the
interior of both the bars, `L1` and `L2`. We will use `interior_bc` method to sample points in the interior of the geometry. Again, the equations to solve are specified as a dictionary input to `outvar_sympy` parameter. These dictionaries are then used when unrolling the computational graph for training.

For this problem we have the `'diffusion_u_1':0` and `'diffusion_u_2':0` for bars `L1` and `L2` respectively. The parameter `bounds`, determines the range for sampling the values for variables $x$ and $y$. The `lambda_sympy` parameter is used to determine the weights for different losses. In this problem, we weight each point equally and hence keep the variable to 1 for each key (default). 


```python
class DiffusionTrain(TrainDomain):
  def __init__(self, **config):
    super(DiffusionTrain, self).__init__()
    # sympy variables
    x = Symbol('x')
    c = Symbol('c')
    
    # right hand side (x = 2) Pt c
    IC = L2.boundary_bc(outvar_sympy={'u_2': Tc},
                        batch_size_per_area=1,
                        criteria=Eq(x, 2))
    self.add(IC, name="RightHandSide")
    
    # left hand side (x = 0) Pt a
    IC = L1.boundary_bc(outvar_sympy={'u_1': Ta},
                        batch_size_per_area=1,
                        criteria=Eq(x, 0))
    self.add(IC, name="LeftHandSide")
    
    # interface 1-2
    IC = L1.boundary_bc(outvar_sympy={'diffusion_interface_dirichlet_u_1_u_2': 0,
                                      'diffusion_interface_neumann_u_1_u_2': 0},
                        lambda_sympy={'lambda_diffusion_interface_dirichlet_u_1_u_2': 1,
                                      'lambda_diffusion_interface_neumann_u_1_u_2': 1},
                        batch_size_per_area=1,
                        criteria=Eq(x, 1))
    self.add(IC, name="Interface1n2")
    
    # interior 1
    interior = L1.interior_bc(outvar_sympy={'diffusion_u_1': 0},
                              lambda_sympy={'lambda_diffusion_u_1': 1},
                              bounds={x: (0, 1)},
                              batch_size_per_area=200)
    self.add(interior, name="Interior1")
    
    # interior 2
    interior = L2.interior_bc(outvar_sympy={'diffusion_u_2': 0},
                              lambda_sympy={'lambda_diffusion_u_2': 1},
                              bounds={x: (1, 2)},
                              batch_size_per_area=200)
    self.add(interior, name="Interior2")
```

At this point you might be wondering where do we input the parameters of the equation, for eg. values for $D_1, D_2$, etc. Don't worry, we will discuss them while making the neural network solver. But before that, let's create the validation data to verify our simulation results against the analytical solution.

## Creating Validation Domain

For this 1d bar problem where the conductivity is constant in each bar, the temperature varies linearly with position inside the solid. The analytical solution can then be given as:

$$
\begin{align}
U_1 = xT_b + (1-x)T_a, && \text{when } 0 \leq x \leq 1 \\
U_2 = (x-1)T_c + (2-x)T_b, && \text{when } 1 \leq x \leq 2 \\
\end{align}
$$

where, 
$$
\begin{align}
T_a = U_1|_{x=0}, && T_c = U_2|_{x=2}, && \frac{\left(T_c + \left( D_1/D_2 \right)T_a \right)}{1+ \left( D_1/D_2 \right)}\\
\end{align}
$$

Now let's create the validation domains. The validation domain is created by inheriting from the `ValidationDomain` parent class. We use numpy to solve for the `u_1` and `u_2` based on the analytical expressions we showed above. The dictionary of generated numpy arrays (`invar_numpy` and `outvar_numpy`)for input and output variables is used as an input to the class method `from_numpy`.


```python
class DiffusionVal(ValidationDomain):
  def __init__(self, **config):
    super(DiffusionVal, self).__init__()

    x = np.expand_dims(np.linspace(0, 1, 100), axis=-1)
    u_1 = x*Tb + (1-x)*Ta
    invar_numpy = {'x': x}
    outvar_numpy = {'u_1': u_1}
    val = Validation.from_numpy(invar_numpy, outvar_numpy)
    self.add(val, name='Val1')
    
    # make validation data line 2
    x = np.expand_dims(np.linspace(1, 2, 100), axis=-1)
    u_2 = (x-1)*Tc + (2-x)*Tb
    invar_numpy = {'x': x}
    outvar_numpy = {'u_2': u_2}
    val = Validation.from_numpy(invar_numpy, outvar_numpy)
    self.add(val, name='Val2')
```

## Creating Monitor Domain

SimNet library allows you to monitor desired quantities in Tensorboard as the simulation progresses and
assess the convergence. A `MonitorDomain` can be used to create such an feature. This a useful feature when we want
to monitor convergence based on a quantity of interest. Examples of such quantities can be point values of variables,
surface averages, volume averages or any other derived quantities. The variables are available as TensorFlow tensors. We can perform tensor operations available in TensorFlow to compute any desired derived quantity of our choice. 

In the code below, we create monitors for flux at the interface. The variable `u_1__x` represents the derivative of `u_1` in x-direction (two underscores (`__`) and the variable (`x`)). The same notation is used while handling other derivatives using the SimNet library. (eg. a neumann boundary condition of $\frac{dU_1}{dx}=0$ can be assigned as `'u_1__x':0` in the train domain for solving the same problem with a adiabatic/fixed flux boundary condition). 

The points to sample can be selected in a similar way as we did for specifying the Train domain. We create the monitors by inheriting from the `MonitorDoamin` parent class


```python
class DiffusionMonitor(MonitorDomain):
  def __init__(self, **config):
    super(DiffusionMonitor, self).__init__()
    x = Symbol('x')

    # flux in U1 at x = 1 
    fluxU1 = Monitor(L1.sample_boundary(10, criteria=Eq(x, 1)),
                    {'flux_U1': lambda var: tf.reduce_mean(D1*var['u_1__x'])})
    self.add(fluxU1, 'FluxU1')

    # flux in U2 at x = 1 
    fluxU2 = Monitor(L2.sample_boundary(10, criteria=Eq(x, 1)),
                    {'flux_U2': lambda var: tf.reduce_mean(D2*var['u_2__x'])})
    self.add(fluxU2, 'FluxU2')
```

## Creating the Neural Network Solver

Now that we have the train domain and other validation and monitor domains defined, we can prepare the neural network solver and run the problem. The solver is defined by inheriting the `Solver` parent class. The `train_domain`, `val_domain`, and `monitor_domains` are assigned. The equations to be solved are specified under `self.equations`. Here, we will call the `Diffusion` and `DiffusionInterface` classes we defined earlier to include the PDEs of the problem. Now we will pass the appropriate values for the parameters like the variable name (eg. `T='u_1'`) and also specify the dimensions of the problem (1d and steady).

The inputs and the outputs of the neural network are specified and the nodes of the architecture are made. The default
network architecture is a simple fully connected multi-layer perceptron architecture with *swish* activation function. The network consists of 6 hidden layers with 512 nodes in each layer. Here we are using two separate neural networks for each variable (`u_1` and `u_2`). All these values can be modified through `update_defaults` function. Also, the different architectures in SimNet library can be used (eg. Fourier Net architecture, Radial Basis Neural Network architecture, etc. More details can be found in *SimNet User Guide*). We use the default exponential learning rate decay, set the start learning rate and decay steps and also assign the `'max_steps'` to 5000. 


```python
# Define neural network
class DiffusionSolver(Solver):
  train_domain = DiffusionTrain
  val_domain = DiffusionVal
  monitor_domain = DiffusionMonitor

  def __init__(self, **config):
    super(DiffusionSolver, self).__init__(**config)

    self.equations = (Diffusion(T='u_1', D=D1, dim=1, time=False).make_node()
                      + Diffusion(T='u_2', D=D2, dim=1, time=False).make_node()
                      + DiffusionInterface('u_1', 'u_2', D1, D2, dim=1, time=False).make_node())
    diff_net_u_1 = self.arch.make_node(name='diff_net_u_1',
                                   inputs=['x'],
                                   outputs=['u_1'])
    diff_net_u_2 = self.arch.make_node(name='diff_net_u_2',
                                   inputs=['x'],
                                   outputs=['u_2'])
    self.nets = [diff_net_u_1, diff_net_u_2]

  @classmethod 
  def update_defaults(cls, defaults):
    defaults.update({
        'network_dir': './network_checkpoint_diff',
        'max_steps': 5000,
        'decay_steps': 100,
        'start_lr': 1e-4,
        #'end_lr': 1e-6,
        })
```

Awesome! We have just completed the file set up for the problem using the SimNet library. We are now ready to solve the PDEs using Neural Networks!

Before we can start training, we can make use of Tensorboard for visualizing the loss values and convergence of several other monitors we just created. This can be done inside the jupyter framework by selecting the directory in which the checkpoint will be stored by clicking on the small checkbox next to it. The option to launch a Tensorboard then shows up in that directory. 



Also, SimNet is desinged such that it can accept command line arguments. This causes issues when the code is directly executed through the jupyter notebook. So as a workaround, we will save the code in form of a python script and execute that script inside the jupyter cell. This example is already saved for you and the code block below executes that script `diffusion_bar.py`. You are encouraged to open the script in a different window and go through the code once before executing. Also, feel free to edit the parameters of the model and see its effect on the results.


```python
import os
import sys

sys.path.append('../../source_code/diffusion_1d')

!python ../../source_code/diffusion_1d/diffusion_bar.py
```

## Visualizing the solution

SimNet saves the data in .vtu and .npz format by default. The .npz arrays can be plotted to visualize the output of the simulation. The .npz files that are created are found in the `network_checkpoint*` directory. 

Now let's plot the temperature along the bar for the analytical and the neural network solution. A sample script to plot the results is shown below. If the training is complete, you should get the results like shown below. As we can see, our neural network solution and the analytical solution match almost exactly for this diffusion problem. 




```python
%%capture
import sys
!{sys.executable} -m pip install ipympl
%matplotlib inline
import matplotlib.pyplot as plt

plt.figure()
network_dir = './network_checkpoint_diff/val_domain/results/'
u_1_pred = np.load(network_dir + 'Val1_pred.npz', allow_pickle=True)
u_2_pred = np.load(network_dir + 'Val2_pred.npz', allow_pickle=True)
u_1_pred = np.atleast_1d(u_1_pred.f.arr_0)[0]
u_2_pred = np.atleast_1d(u_2_pred.f.arr_0)[0]

plt.plot(u_1_pred['x'][:,0], u_1_pred['u_1'][:,0], '--', label='u_1_pred')
plt.plot(u_2_pred['x'][:,0], u_2_pred['u_2'][:,0], '--', label='u_2_pred')

u_1_true = np.load(network_dir + 'Val1_true.npz', allow_pickle=True)
u_2_true = np.load(network_dir + 'Val2_true.npz', allow_pickle=True)
u_1_true = np.atleast_1d(u_1_true.f.arr_0)[0]
u_2_true = np.atleast_1d(u_2_true.f.arr_0)[0]

plt.plot(u_1_true['x'][:,0], u_1_true['u_1'][:,0], label='u_1_true')
plt.plot(u_2_true['x'][:,0], u_2_true['u_2'][:,0], label='u_2_true')

plt.legend()

plt.savefig('image_diffusion_problem_bootcamp.png')
```


```python
from IPython.display import Image
Image(filename='image_diffusion_problem_bootcamp.png') 
```



# Parameterizing the PDE

As we discussed in the introductory notebook, one important advantage of a PINN solver over traditional numerical methods is its ability to solve parameterized geometries and PDEs. This was initially proposed in the [paper](https://arxiv.org/abs/1906.02382) published by Sun et al. This allows us to significant computational advantage as one can now use PINNs to solve for multiple desigs/cases in a single training. Once the training is complete, it is possible to run inference on several geometry/physical parameter combinations as a post-processing step, without solving the forward problem again. 

To demonstrate the concept, we will train the same 1d diffusion problem, but now by parameterizing the conductivity of the first bar in the range (5, 25). Once the training is complete, we can obtain the results for any conductivity value in that range saving us the time to train multiple models. 

## Case Setup

The definition of equations remain the same for this part. Since earlier while defining the equations, we already defined the constants and coefficients of the PDE to be either numerical values or strings, this will allow us to paramterize `D1` by passing it as a string while calling the equations and making the neural network. Now let's start by creating the paramterized train domain. We will skip the parts that are common to the previous section and only discuss the changes. The complete script can be referred in `diffusion_bar_parameterized.py`

## Creating Train Domain

Before starting out to create the train domain using the `TrainDomain` class, we will create the symbolic variable for the $D_1$ and also specify the range of variation for the variable. While the simulation runs, we will validate it against the same diffusion coefficient that we solved earlier i.e. $D_1=10$. Once the sybolic variables and the ranges are described for sampling, these parameter ranges need to be inputted to the `param_ranges` attribute of each boundary and internal sampling (`boundary_bc` and `interior_bc`)


```python
# params for domain
L1 = Line1D(0,1)
L2 = Line1D(1,2)

D1 = Symbol('D1')
D1_range = {D1: (5, 25)}
D1_validation = 1e1

D2 = 1e-1

Tc = 100
Ta = 0
Tb = (Tc + (D1/D2)*Ta)/(1 + (D1/D2))

Tb_validation = float(Tb.evalf(subs={D1: 1e1}))

class DiffusionTrain(TrainDomain):
  def __init__(self, **config):
    super(DiffusionTrain, self).__init__()
    # sympy variables
    x = Symbol('x')
    c = Symbol('c')

    # right hand side (x = 2) Pt c
    IC = L2.boundary_bc(outvar_sympy={'u_2': Tc},
                        batch_size_per_area=10,
                        criteria=Eq(x, 2),
                        param_ranges=D1_range)
    self.add(IC, name="RightHandSide")

    # left hand side (x = 0) Pt a
    IC = L1.boundary_bc(outvar_sympy={'u_1': Ta},
                        batch_size_per_area=10,
                        criteria=Eq(x, 0),
                        param_ranges=D1_range)
    self.add(IC, name="LeftHandSide")

    # interface 1-2
    IC = L1.boundary_bc(outvar_sympy={'diffusion_interface_dirichlet_u_1_u_2': 0,
                                      'diffusion_interface_neumann_u_1_u_2': 0},
                        lambda_sympy={'lambda_diffusion_interface_dirichlet_u_1_u_2': 1,
                                      'lambda_diffusion_interface_neumann_u_1_u_2': 1},
                        batch_size_per_area=10,
                        criteria=Eq(x, 1),
                        param_ranges=D1_range)
    self.add(IC, name="Interface1n2")

    # interior 1
    interior = L1.interior_bc(outvar_sympy={'diffusion_u_1': 0},
                              lambda_sympy={'lambda_diffusion_u_1': 1},
                              bounds={x: (0, 1)},
                              batch_size_per_area=400,
                              param_ranges=D1_range)
    self.add(interior, name="Interior1")

    # interior 2
    interior = L2.interior_bc(outvar_sympy={'diffusion_u_2': 0},
                              lambda_sympy={'lambda_diffusion_u_2': 1},
                              bounds={x: (1, 2)},
                              batch_size_per_area=400,
                              param_ranges=D1_range)
    self.add(interior, name="Interior2")
```

## Creating Validation and Monitor Domains

The process to create these domains is again similar to the previous section. For validation data, we need to create an additional key for the string `'D1'` in the `invar_numpy`. The value for this string can be in the range we specified earlier and which we would like to validate against. It is possible to create multiple validations if required, eg. different $D_1$ values. For monitor domain, similar to `interior_bc` and `boundary_bc` in the train domain, we will supply the paramter ranges for monitoring in the `param_ranges` attribute of the `sample_boundary` method. 


```python
class DiffusionVal(ValidationDomain):
  def __init__(self, **config):
    super(DiffusionVal, self).__init__()
    # make validation data line 1
    x = np.expand_dims(np.linspace(0, 1, 100), axis=-1)
    D1 = np.zeros_like(x) + D1_validation                      # For creating D1 input array
    u_1 = x*Tb_validation + (1-x)*Ta
    invar_numpy = {'x': x}                                    # Set the invars for the required D1
    invar_numpy.update({'D1': np.full_like(invar_numpy['x'], D1_validation)})
    outvar_numpy = {'u_1': u_1}
    val = Validation.from_numpy(invar_numpy, outvar_numpy)
    self.add(val, name='Val1')

    # make validation data line 2
    x = np.expand_dims(np.linspace(1, 2, 100), axis=-1)
    u_2 = (x-1)*Tc + (2-x)*Tb_validation
    invar_numpy = {'x': x}                           # Set the invars for the required D1
    invar_numpy.update({'D1': np.full_like(invar_numpy['x'], D1_validation)})
    outvar_numpy = {'u_2': u_2}
    val = Validation.from_numpy(invar_numpy, outvar_numpy)
    self.add(val, name='Val2')

class DiffusionMonitor(MonitorDomain):
  def __init__(self, **config):
    super(DiffusionMonitor, self).__init__()
    x = Symbol('x')

    # flux in U1 at x = 1
    fluxU1 = Monitor(L1.sample_boundary(10, criteria=Eq(x, 1), param_ranges={D1: D1_validation}),   # Set the parameter range for the required D1
                    {'flux_U1': lambda var: tf.reduce_mean(D1_validation*var['u_1__x'])})
    self.add(fluxU1, 'FluxU1')

    # flux in U2 at x = 1
    fluxU2 = Monitor(L2.sample_boundary(10, criteria=Eq(x, 1), param_ranges={D1: D1_validation}),   # Set the parameter range for the required D1
                    {'flux_U2': lambda var: tf.reduce_mean(D2*var['u_2__x'])})
    self.add(fluxU2, 'FluxU2')
```

## Creating the Neural Network Solver

Once all the parameterized domain definitions are completed, for training the parameterized model, we will have the symbolic parameters we defined earlier as inputs to both the neural networks in `diff_net_u_1` and `diff_net_u_2` viz. `'D1'` along with the usual x coordinate. The outputs remain the same as what we would have for any other non-parameterized simulation. 


```python
# Define neural network
class DiffusionSolver(Solver):
  train_domain = DiffusionTrain
  val_domain = DiffusionVal
  monitor_domain = DiffusionMonitor

  def __init__(self, **config):
    super(DiffusionSolver, self).__init__(**config)

    self.equations = (Diffusion(T='u_1', D='D1', dim=1, time=False).make_node()         # Symbolic input to the equation
                      + Diffusion(T='u_2', D=D2, dim=1, time=False).make_node()
                      + DiffusionInterface('u_1', 'u_2', 'D1', D2, dim=1, time=False).make_node())
    diff_net_u_1 = self.arch.make_node(name='diff_net_u_1',
                                   inputs=['x', 'D1'],                                  # Add the parameters to the network
                                   outputs=['u_1'])
    diff_net_u_2 = self.arch.make_node(name='diff_net_u_2',
                                   inputs=['x', 'D1'],
                                   outputs=['u_2'])
    self.nets = [diff_net_u_1, diff_net_u_2]

  @classmethod # Explain This
  def update_defaults(cls, defaults):
    defaults.update({
        'network_dir': './network_checkpoint_diff_parameterized',
        'max_steps': 10000,
        'decay_steps': 200,
        'start_lr': 1e-4,
        'layer_size': 256,
        'xla': True,
        })
```

## Visualizing the solution

The .npz arrays can be plotted similar to previous section to visualize the output of the simulation. You can see that we get the same answer as the analytical solution. You can try to run the problem in `eval` mode by chaning the validation data and see how it performs for the other `D1` values as well. To run the model in evaluation mode (i.e. without training), you just need to add the `--run_mode=eval` flag while executing the script. 



You can see that at a fractional increase in computational time, we solved the PDE for $D_1$ ranging from (5, 25). This concept can easily be extended to more complicated problems and this ability of solving parameterized problems comes very handy during desing optimization and exploring the desing space. For more examples of solving parameterized problems, please refer to *SimNet User Guide Chapter 13*

# Licensing
This material is released by NVIDIA Corporation under the Creative Commons Attribution 4.0 International (CC BY 4.0)

&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&ensp;
[Home Page](../../Start_Here.ipynb)

[Previous Notebook](../introduction/Introductory_Notebook.ipynb)
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&ensp;
[1](../introduction/Introductory_Notebook.ipynb)
[2]
[3](../spring_mass/Spring_Mass_Problem_Notebook.ipynb)
[4](../chip_2d/Challenge_CFD_Problem_Notebook.ipynb)
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
&emsp;&emsp;&emsp;&emsp;&emsp;
[Next Notebook](../spring_mass/Spring_Mass_Problem_Notebook.ipynb)

42932f2a831b29796a47b33149b95317a022112c
ipynb
Jupyter Notebook
hpc_ai/PINN/English/python/jupyter_notebook/diffusion_1d/Diffusion_Problem_Notebook.ipynb
AditiRM/gpubootcamp
bd14fec2383d7f116760ca9f4efd5816523f9be9
Apache-2.0
2020-11-19T21:12:48.000Z
2022-03-31T13:47:01.000Z
hpc_ai/PINN/English/python/jupyter_notebook/diffusion_1d/Diffusion_Problem_Notebook.ipynb
bharatk-parallel/gpubootcamp-1
7cce9760fdc6f211790ce1d320234acda0aa8808
Apache-2.0
2021-05-05T17:03:18.000Z
2022-03-31T07:23:36.000Z
hpc_ai/PINN/English/python/jupyter_notebook/diffusion_1d/Diffusion_Problem_Notebook.ipynb
bharatk-parallel/gpubootcamp-1
7cce9760fdc6f211790ce1d320234acda0aa8808
Apache-2.0
2020-11-19T21:32:21.000Z
2022-03-31T11:08:41.000Z
Qwen/Qwen-72B
1. YES
2. YES