update examples

This commit is contained in:
Dario Coscia
2023-05-29 15:34:23 +02:00
committed by Nicola Demo
parent 37e9658211
commit eb531747e5
11 changed files with 331 additions and 216 deletions

View File

@@ -5,13 +5,26 @@ from pina.operators import grad
from pina import Condition from pina import Condition
from pina.span import Span from pina.span import Span
# ===================================================== #
# #
# This script implements the one dimensional Burger #
# problem. The Burgers1D class is defined inheriting #
# from TimeDependentProblem, SpatialProblem and we #
# denote: #
# u --> field variable #
# x --> spatial variable #
# t --> temporal variable #
# #
# ===================================================== #
class Burgers1D(TimeDependentProblem, SpatialProblem): class Burgers1D(TimeDependentProblem, SpatialProblem):
# assign output/ spatial and temporal variables
output_variables = ['u'] output_variables = ['u']
spatial_domain = Span({'x': [-1, 1]}) spatial_domain = Span({'x': [-1, 1]})
temporal_domain = Span({'t': [0, 1]}) temporal_domain = Span({'t': [0, 1]})
# define the burger equation
def burger_equation(input_, output_): def burger_equation(input_, output_):
du = grad(output_, input_) du = grad(output_, input_)
ddu = grad(du, input_, components=['dudx']) ddu = grad(du, input_, components=['dudx'])
@@ -21,17 +34,20 @@ class Burgers1D(TimeDependentProblem, SpatialProblem):
(0.01/torch.pi)*ddu.extract(['ddudxdx']) (0.01/torch.pi)*ddu.extract(['ddudxdx'])
) )
# define nill dirichlet boundary conditions
def nil_dirichlet(input_, output_): def nil_dirichlet(input_, output_):
u_expected = 0.0 u_expected = 0.0
return output_.extract(['u']) - u_expected return output_.extract(['u']) - u_expected
# define initial condition
def initial_condition(input_, output_): def initial_condition(input_, output_):
u_expected = -torch.sin(torch.pi*input_.extract(['x'])) u_expected = -torch.sin(torch.pi*input_.extract(['x']))
return output_.extract(['u']) - u_expected return output_.extract(['u']) - u_expected
# problem condition statement
conditions = { conditions = {
'gamma1': Condition(Span({'x': -1, 't': [0, 1]}), nil_dirichlet), 'gamma1': Condition(location=Span({'x': -1, 't': [0, 1]}), function=nil_dirichlet),
'gamma2': Condition(Span({'x': 1, 't': [0, 1]}), nil_dirichlet), 'gamma2': Condition(location=Span({'x': 1, 't': [0, 1]}), function=nil_dirichlet),
't0': Condition(Span({'x': [-1, 1], 't': 0}), initial_condition), 't0': Condition(location=Span({'x': [-1, 1], 't': 0}), function=initial_condition),
'D': Condition(Span({'x': [-1, 1], 't': [0, 1]}), burger_equation), 'D': Condition(location=Span({'x': [-1, 1], 't': [0, 1]}), function=burger_equation),
} }

View File

@@ -1,45 +1,47 @@
import torch # import torch
from pina.problem import Problem # from pina.problem import Problem
from pina.segment import Segment # from pina.segment import Segment
from pina.cube import Cube # from pina.cube import Cube
from pina.problem2d import Problem2D # from pina.problem2d import Problem2D
xmin, xmax, ymin, ymax = -1, 1, -1, 1 # xmin, xmax, ymin, ymax = -1, 1, -1, 1
class EllipticOptimalControl(Problem2D): # class EllipticOptimalControl(Problem2D):
def __init__(self, alpha=1): # def __init__(self, alpha=1):
def term1(input_, output_): # def term1(input_, output_):
grad_p = self.grad(output_.extract(['p']), input_) # grad_p = self.grad(output_.extract(['p']), input_)
gradgrad_p_x1 = self.grad(grad_p.extract(['x1']), input_) # gradgrad_p_x1 = self.grad(grad_p.extract(['x1']), input_)
gradgrad_p_x2 = self.grad(grad_p.extract(['x2']), input_) # gradgrad_p_x2 = self.grad(grad_p.extract(['x2']), input_)
yd = 2.0 # yd = 2.0
return output_.extract(['y']) - yd - (gradgrad_p_x1.extract(['x1']) + gradgrad_p_x2.extract(['x2'])) # return output_.extract(['y']) - yd - (gradgrad_p_x1.extract(['x1']) + gradgrad_p_x2.extract(['x2']))
def term2(input_, output_): # def term2(input_, output_):
grad_y = self.grad(output_.extract(['y']), input_) # grad_y = self.grad(output_.extract(['y']), input_)
gradgrad_y_x1 = self.grad(grad_y.extract(['x1']), input_) # gradgrad_y_x1 = self.grad(grad_y.extract(['x1']), input_)
gradgrad_y_x2 = self.grad(grad_y.extract(['x2']), input_) # gradgrad_y_x2 = self.grad(grad_y.extract(['x2']), input_)
return - (gradgrad_y_x1.extract(['x1']) + gradgrad_y_x2.extract(['x2'])) - output_.extract(['u']) # return - (gradgrad_y_x1.extract(['x1']) + gradgrad_y_x2.extract(['x2'])) - output_.extract(['u'])
def term3(input_, output_): # def term3(input_, output_):
return output_.extract(['p']) - output_.extract(['u'])*alpha # return output_.extract(['p']) - output_.extract(['u'])*alpha
def nil_dirichlet(input_, output_): # def nil_dirichlet(input_, output_):
y_value = 0.0 # y_value = 0.0
p_value = 0.0 # p_value = 0.0
return torch.abs(output_.extract(['y']) - y_value) + torch.abs(output_.extract(['p']) - p_value) # return torch.abs(output_.extract(['y']) - y_value) + torch.abs(output_.extract(['p']) - p_value)
self.conditions = { # self.conditions = {
'gamma1': {'location': Segment((xmin, ymin), (xmax, ymin)), 'func': nil_dirichlet}, # 'gamma1': {'location': Segment((xmin, ymin), (xmax, ymin)), 'func': nil_dirichlet},
'gamma2': {'location': Segment((xmax, ymin), (xmax, ymax)), 'func': nil_dirichlet}, # 'gamma2': {'location': Segment((xmax, ymin), (xmax, ymax)), 'func': nil_dirichlet},
'gamma3': {'location': Segment((xmax, ymax), (xmin, ymax)), 'func': nil_dirichlet}, # 'gamma3': {'location': Segment((xmax, ymax), (xmin, ymax)), 'func': nil_dirichlet},
'gamma4': {'location': Segment((xmin, ymax), (xmin, ymin)), 'func': nil_dirichlet}, # 'gamma4': {'location': Segment((xmin, ymax), (xmin, ymin)), 'func': nil_dirichlet},
'D1': {'location': Cube([[xmin, xmax], [ymin, ymax]]), 'func': [term1, term2, term3]}, # 'D1': {'location': Cube([[xmin, xmax], [ymin, ymax]]), 'func': [term1, term2, term3]},
} # }
self.input_variables = ['x1', 'x2'] # self.input_variables = ['x1', 'x2']
self.output_variables = ['u', 'p', 'y'] # self.output_variables = ['u', 'p', 'y']
self.spatial_domain = Cube([[xmin, xmax], [xmin, xmax]]) # self.spatial_domain = Cube([[xmin, xmax], [xmin, xmax]])
raise NotImplementedError('not available problem at the moment...')

View File

@@ -0,0 +1,46 @@
from pina.problem import SpatialProblem
from pina import Condition, Span
from pina.operators import grad
import torch
# ===================================================== #
# #
# This script implements a simple first order ode. #
# The FirstOrderODE class is defined inheriting from #
# SpatialProblem. We denote: #
# y --> field variable #
# x --> spatial variable #
# #
# ===================================================== #
class FirstOrderODE(SpatialProblem):
# variable domain range
x_rng = [0, 5]
# field variable
output_variables = ['y']
# create domain
spatial_domain = Span({'x': x_rng})
# define the ode
def ode(input_, output_):
y = output_
x = input_
return grad(y, x) + y - x
# define initial conditions
def fixed(input_, output_):
exp_value = 1.
return output_ - exp_value
# define real solution
def solution(self, input_):
x = input_
return x - 1.0 + 2*torch.exp(-x)
# define problem conditions
conditions = {
'bc': Condition(location=Span({'x': x_rng[0]}), function=fixed),
'dd': Condition(location=Span({'x': x_rng}), function=ode),
}
truth_solution = solution

View File

@@ -1,53 +1,53 @@
import numpy as np # import numpy as np
import torch # import torch
from pina.problem import Problem # from pina.problem import Problem
from pina.segment import Segment # from pina.segment import Segment
from pina.cube import Cube # from pina.cube import Cube
from pina.problem2d import Problem2D # from pina.problem2d import Problem2D
xmin, xmax, ymin, ymax = -1, 1, -1, 1 # xmin, xmax, ymin, ymax = -1, 1, -1, 1
class ParametricEllipticOptimalControl(Problem2D): # class ParametricEllipticOptimalControl(Problem2D):
def __init__(self, alpha=1): # def __init__(self, alpha=1):
def term1(input_, param_, output_): # def term1(input_, param_, output_):
grad_p = self.grad(output_['p'], input_) # grad_p = self.grad(output_['p'], input_)
gradgrad_p_x1 = self.grad(grad_p['x1'], input_) # gradgrad_p_x1 = self.grad(grad_p['x1'], input_)
gradgrad_p_x2 = self.grad(grad_p['x2'], input_) # gradgrad_p_x2 = self.grad(grad_p['x2'], input_)
return output_['y'] - param_ - (gradgrad_p_x1['x1'] + gradgrad_p_x2['x2']) # return output_['y'] - param_ - (gradgrad_p_x1['x1'] + gradgrad_p_x2['x2'])
def term2(input_, param_, output_): # def term2(input_, param_, output_):
grad_y = self.grad(output_['y'], input_) # grad_y = self.grad(output_['y'], input_)
gradgrad_y_x1 = self.grad(grad_y['x1'], input_) # gradgrad_y_x1 = self.grad(grad_y['x1'], input_)
gradgrad_y_x2 = self.grad(grad_y['x2'], input_) # gradgrad_y_x2 = self.grad(grad_y['x2'], input_)
return - (gradgrad_y_x1['x1'] + gradgrad_y_x2['x2']) - output_['u_param'] # return - (gradgrad_y_x1['x1'] + gradgrad_y_x2['x2']) - output_['u_param']
def term3(input_, param_, output_): # def term3(input_, param_, output_):
return output_['p'] - output_['u_param']*alpha # return output_['p'] - output_['u_param']*alpha
def term(input_, param_, output_): # def term(input_, param_, output_):
return term1( input_, param_, output_) +term2( input_, param_, output_) + term3( input_, param_, output_) # return term1( input_, param_, output_) +term2( input_, param_, output_) + term3( input_, param_, output_)
def nil_dirichlet(input_, param_, output_): # def nil_dirichlet(input_, param_, output_):
y_value = 0.0 # y_value = 0.0
p_value = 0.0 # p_value = 0.0
return torch.abs(output_['y'] - y_value) + torch.abs(output_['p'] - p_value) # return torch.abs(output_['y'] - y_value) + torch.abs(output_['p'] - p_value)
self.conditions = { # self.conditions = {
'gamma1': {'location': Segment((xmin, ymin), (xmax, ymin)), 'func': nil_dirichlet}, # 'gamma1': {'location': Segment((xmin, ymin), (xmax, ymin)), 'func': nil_dirichlet},
'gamma2': {'location': Segment((xmax, ymin), (xmax, ymax)), 'func': nil_dirichlet}, # 'gamma2': {'location': Segment((xmax, ymin), (xmax, ymax)), 'func': nil_dirichlet},
'gamma3': {'location': Segment((xmax, ymax), (xmin, ymax)), 'func': nil_dirichlet}, # 'gamma3': {'location': Segment((xmax, ymax), (xmin, ymax)), 'func': nil_dirichlet},
'gamma4': {'location': Segment((xmin, ymax), (xmin, ymin)), 'func': nil_dirichlet}, # 'gamma4': {'location': Segment((xmin, ymax), (xmin, ymin)), 'func': nil_dirichlet},
'D1': {'location': Cube([[xmin, xmax], [ymin, ymax]]), 'func': term}, # 'D1': {'location': Cube([[xmin, xmax], [ymin, ymax]]), 'func': term},
#'D2': {'location': Cube([[0, 1], [0, 1]]), 'func': term2}, # #'D2': {'location': Cube([[0, 1], [0, 1]]), 'func': term2},
#'D3': {'location': Cube([[0, 1], [0, 1]]), 'func': term3} # #'D3': {'location': Cube([[0, 1], [0, 1]]), 'func': term3}
} # }
self.input_variables = ['x1', 'x2']
self.output_variables = ['u', 'p', 'y']
self.parameters = ['mu']
self.spatial_domain = Cube([[xmin, xmax], [xmin, xmax]])
self.parameter_domain = np.array([[0.5, 3]])
# self.input_variables = ['x1', 'x2']
# self.output_variables = ['u', 'p', 'y']
# self.parameters = ['mu']
# self.spatial_domain = Cube([[xmin, xmax], [xmin, xmax]])
# self.parameter_domain = np.array([[0.5, 3]])
raise NotImplementedError('not available problem at the moment...')

View File

@@ -5,22 +5,42 @@ from pina import Span, Condition
from pina.problem import SpatialProblem, ParametricProblem from pina.problem import SpatialProblem, ParametricProblem
from pina.operators import grad, nabla from pina.operators import grad, nabla
# ===================================================== #
# #
# This script implements the two dimensional #
# Parametric Elliptic Optimal Control problem. #
# The ParametricEllipticOptimalControl class is #
# inherited from TimeDependentProblem, SpatialProblem #
# and we denote: #
# u --> field variable #
# p --> field variable #
# y --> field variable #
# x1, x2 --> spatial variables #
# mu, alpha --> problem parameters #
# #
# More info in https://arxiv.org/pdf/2110.13530.pdf #
# Section 4.2 of the article #
# ===================================================== #
class ParametricEllipticOptimalControl(SpatialProblem, ParametricProblem): class ParametricEllipticOptimalControl(SpatialProblem, ParametricProblem):
# setting spatial variables ranges
xmin, xmax, ymin, ymax = -1, 1, -1, 1 xmin, xmax, ymin, ymax = -1, 1, -1, 1
x_range = [xmin, xmax]
y_range = [ymin, ymax]
# setting parameters range
amin, amax = 0.0001, 1 amin, amax = 0.0001, 1
mumin, mumax = 0.5, 3 mumin, mumax = 0.5, 3
mu_range = [mumin, mumax] mu_range = [mumin, mumax]
a_range = [amin, amax] a_range = [amin, amax]
x_range = [xmin, xmax] # setting field variables
y_range = [ymin, ymax]
output_variables = ['u', 'p', 'y'] output_variables = ['u', 'p', 'y']
# setting spatial and parameter domain
spatial_domain = Span({'x1': x_range, 'x2': y_range}) spatial_domain = Span({'x1': x_range, 'x2': y_range})
parameter_domain = Span({'mu': mu_range, 'alpha': a_range}) parameter_domain = Span({'mu': mu_range, 'alpha': a_range})
# equation terms as in https://arxiv.org/pdf/2110.13530.pdf
def term1(input_, output_): def term1(input_, output_):
laplace_p = nabla(output_, input_, components=['p'], d=['x1', 'x2']) laplace_p = nabla(output_, input_, components=['p'], d=['x1', 'x2'])
return output_.extract(['y']) - input_.extract(['mu']) - laplace_p return output_.extract(['y']) - input_.extract(['mu']) - laplace_p
@@ -37,21 +57,22 @@ class ParametricEllipticOptimalControl(SpatialProblem, ParametricProblem):
p_exp = 0.0 p_exp = 0.0
return output_.extract(['p']) - p_exp return output_.extract(['p']) - p_exp
# setting problem condition formulation
conditions = { conditions = {
'gamma1': Condition( 'gamma1': Condition(
Span({'x1': x_range, 'x2': 1, 'mu': mu_range, 'alpha': a_range}), location=Span({'x1': x_range, 'x2': 1, 'mu': mu_range, 'alpha': a_range}),
[state_dirichlet, adj_dirichlet]), function=[state_dirichlet, adj_dirichlet]),
'gamma2': Condition( 'gamma2': Condition(
Span({'x1': x_range, 'x2': -1, 'mu': mu_range, 'alpha': a_range}), location=Span({'x1': x_range, 'x2': -1, 'mu': mu_range, 'alpha': a_range}),
[state_dirichlet, adj_dirichlet]), function=[state_dirichlet, adj_dirichlet]),
'gamma3': Condition( 'gamma3': Condition(
Span({'x1': 1, 'x2': y_range, 'mu': mu_range, 'alpha': a_range}), location=Span({'x1': 1, 'x2': y_range, 'mu': mu_range, 'alpha': a_range}),
[state_dirichlet, adj_dirichlet]), function=[state_dirichlet, adj_dirichlet]),
'gamma4': Condition( 'gamma4': Condition(
Span({'x1': -1, 'x2': y_range, 'mu': mu_range, 'alpha': a_range}), location=Span({'x1': -1, 'x2': y_range, 'mu': mu_range, 'alpha': a_range}),
[state_dirichlet, adj_dirichlet]), function=[state_dirichlet, adj_dirichlet]),
'D': Condition( 'D': Condition(
Span({'x1': x_range, 'x2': y_range, location=Span({'x1': x_range, 'x2': y_range,
'mu': mu_range, 'alpha': a_range}), 'mu': mu_range, 'alpha': a_range}),
[term1, term2]), function=[term1, term2]),
} }

View File

@@ -4,37 +4,52 @@ from pina.problem import SpatialProblem, ParametricProblem
from pina.operators import nabla from pina.operators import nabla
from pina import Span, Condition from pina import Span, Condition
# ===================================================== #
# #
# This script implements the two dimensional #
# Parametric Poisson problem. The ParametricPoisson #
# class is defined inheriting from SpatialProblem and #
# ParametricProblem. We denote: #
# u --> field variable #
# x,y --> spatial variables #
# mu1, mu2 --> parameter variables #
# #
# ===================================================== #
class ParametricPoisson(SpatialProblem, ParametricProblem): class ParametricPoisson(SpatialProblem, ParametricProblem):
# assign output/ spatial and parameter variables
output_variables = ['u'] output_variables = ['u']
spatial_domain = Span({'x': [-1, 1], 'y': [-1, 1]}) spatial_domain = Span({'x': [-1, 1], 'y': [-1, 1]})
parameter_domain = Span({'mu1': [-1, 1], 'mu2': [-1, 1]}) parameter_domain = Span({'mu1': [-1, 1], 'mu2': [-1, 1]})
# define the laplace equation
def laplace_equation(input_, output_): def laplace_equation(input_, output_):
force_term = torch.exp( force_term = torch.exp(
- 2*(input_.extract(['x']) - input_.extract(['mu1']))**2 - 2*(input_.extract(['x']) - input_.extract(['mu1']))**2
- 2*(input_.extract(['y']) - input_.extract(['mu2']))**2) - 2*(input_.extract(['y']) - input_.extract(['mu2']))**2)
return nabla(output_.extract(['u']), input_) - force_term return nabla(output_.extract(['u']), input_) - force_term
# define nill dirichlet boundary conditions
def nil_dirichlet(input_, output_): def nil_dirichlet(input_, output_):
value = 0.0 value = 0.0
return output_.extract(['u']) - value return output_.extract(['u']) - value
# problem condition statement
conditions = { conditions = {
'gamma1': Condition( 'gamma1': Condition(
Span({'x': [-1, 1], 'y': 1, 'mu1': [-1, 1], 'mu2': [-1, 1]}), location=Span({'x': [-1, 1], 'y': 1, 'mu1': [-1, 1], 'mu2': [-1, 1]}),
nil_dirichlet), function=nil_dirichlet),
'gamma2': Condition( 'gamma2': Condition(
Span({'x': [-1, 1], 'y': -1, 'mu1': [-1, 1], 'mu2': [-1, 1]}), location=Span({'x': [-1, 1], 'y': -1, 'mu1': [-1, 1], 'mu2': [-1, 1]}),
nil_dirichlet), function=nil_dirichlet),
'gamma3': Condition( 'gamma3': Condition(
Span({'x': 1, 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}), location=Span({'x': 1, 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}),
nil_dirichlet), function=nil_dirichlet),
'gamma4': Condition( 'gamma4': Condition(
Span({'x': -1, 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}), location=Span({'x': -1, 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}),
nil_dirichlet), function=nil_dirichlet),
'D': Condition( 'D': Condition(
Span({'x': [-1, 1], 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}), location=Span({'x': [-1, 1], 'y': [-1, 1], 'mu1': [-1, 1], 'mu2': [-1, 1]}),
laplace_equation), function=laplace_equation),
} }

View File

@@ -5,35 +5,49 @@ from pina.problem import SpatialProblem
from pina.operators import nabla from pina.operators import nabla
from pina import Condition, Span from pina import Condition, Span
# ===================================================== #
# #
# This script implements the two dimensional #
# Poisson problem. The Poisson class is defined #
# inheriting from SpatialProblem. We denote: #
# u --> field variable #
# x,y --> spatial variables #
# #
# ===================================================== #
class Poisson(SpatialProblem): class Poisson(SpatialProblem):
# assign output/ spatial variables
output_variables = ['u'] output_variables = ['u']
spatial_domain = Span({'x': [0, 1], 'y': [0, 1]}) spatial_domain = Span({'x': [0, 1], 'y': [0, 1]})
# define the laplace equation
def laplace_equation(input_, output_): def laplace_equation(input_, output_):
force_term = (torch.sin(input_.extract(['x'])*torch.pi) * force_term = (torch.sin(input_.extract(['x'])*torch.pi) *
torch.sin(input_.extract(['y'])*torch.pi)) torch.sin(input_.extract(['y'])*torch.pi))
nabla_u = nabla(output_.extract(['u']), input_) nabla_u = nabla(output_.extract(['u']), input_)
return nabla_u - force_term return nabla_u - force_term
# define nill dirichlet boundary conditions
def nil_dirichlet(input_, output_): def nil_dirichlet(input_, output_):
value = 0.0 value = 0.0
return output_.extract(['u']) - value return output_.extract(['u']) - value
# problem condition statement
conditions = { conditions = {
'gamma1': Condition(Span({'x': [0, 1], 'y': 1}), nil_dirichlet), 'gamma1': Condition(location=Span({'x': [0, 1], 'y': 1}), function=nil_dirichlet),
'gamma2': Condition(Span({'x': [0, 1], 'y': 0}), nil_dirichlet), 'gamma2': Condition(location=Span({'x': [0, 1], 'y': 0}), function=nil_dirichlet),
'gamma3': Condition(Span({'x': 1, 'y': [0, 1]}), nil_dirichlet), 'gamma3': Condition(location=Span({'x': 1, 'y': [0, 1]}),function=nil_dirichlet),
'gamma4': Condition(Span({'x': 0, 'y': [0, 1]}), nil_dirichlet), 'gamma4': Condition(location=Span({'x': 0, 'y': [0, 1]}), function=nil_dirichlet),
'D': Condition(Span({'x': [0, 1], 'y': [0, 1]}), laplace_equation), 'D': Condition(location=Span({'x': [0, 1], 'y': [0, 1]}), function=laplace_equation),
} }
# real poisson solution
def poisson_sol(self, pts): def poisson_sol(self, pts):
return -( return -(
torch.sin(pts.extract(['x'])*torch.pi)* torch.sin(pts.extract(['x'])*torch.pi)*
torch.sin(pts.extract(['y'])*torch.pi) torch.sin(pts.extract(['y'])*torch.pi)
)/(2*torch.pi**2) )/(2*torch.pi**2)
#return -(np.sin(x*np.pi)*np.sin(y*np.pi))/(2*np.pi**2)
truth_solution = poisson_sol truth_solution = poisson_sol

View File

@@ -5,36 +5,62 @@ from pina.problem import SpatialProblem
from pina.operators import nabla, grad, div from pina.operators import nabla, grad, div
from pina import Condition, Span, LabelTensor from pina import Condition, Span, LabelTensor
# ===================================================== #
# #
# This script implements the two dimensional #
# Stokes problem. The Stokes class is defined #
# inheriting from SpatialProblem. We denote: #
# ux --> field variable velocity along x #
# uy --> field variable velocity along y #
# p --> field variable pressure #
# x,y --> spatial variables #
# #
# ===================================================== #
class Stokes(SpatialProblem): class Stokes(SpatialProblem):
# assign output/ spatial variables
output_variables = ['ux', 'uy', 'p'] output_variables = ['ux', 'uy', 'p']
spatial_domain = Span({'x': [-2, 2], 'y': [-1, 1]}) spatial_domain = Span({'x': [-2, 2], 'y': [-1, 1]})
# define the momentum equation
def momentum(input_, output_): def momentum(input_, output_):
nabla_ = torch.hstack((LabelTensor(nabla(output_.extract(['ux']), input_), ['x']), nabla_ = torch.hstack((LabelTensor(nabla(output_.extract(['ux']), input_), ['x']),
LabelTensor(nabla(output_.extract(['uy']), input_), ['y']))) LabelTensor(nabla(output_.extract(['uy']), input_), ['y'])))
return - nabla_ + grad(output_.extract(['p']), input_) return - nabla_ + grad(output_.extract(['p']), input_)
# define the continuity equation
def continuity(input_, output_): def continuity(input_, output_):
return div(output_.extract(['ux', 'uy']), input_) return div(output_.extract(['ux', 'uy']), input_)
# define the inlet velocity
def inlet(input_, output_): def inlet(input_, output_):
value = 2 * (1 - input_.extract(['y'])**2) value = 2 * (1 - input_.extract(['y'])**2)
return output_.extract(['ux']) - value return output_.extract(['ux']) - value
# define the outlet pressure
def outlet(input_, output_): def outlet(input_, output_):
value = 0.0 value = 0.0
return output_.extract(['p']) - value return output_.extract(['p']) - value
# define the wall condition
def wall(input_, output_): def wall(input_, output_):
value = 0.0 value = 0.0
return output_.extract(['ux', 'uy']) - value return output_.extract(['ux', 'uy']) - value
# problem condition statement
conditions = { conditions = {
'gamma_top': Condition(Span({'x': [-2, 2], 'y': 1}), wall), 'gamma_top': Condition(location=Span({'x': [-2, 2], 'y': 1}), function=wall),
'gamma_bot': Condition(Span({'x': [-2, 2], 'y': -1}), wall), 'gamma_bot': Condition(location=Span({'x': [-2, 2], 'y': -1}), function=wall),
'gamma_out': Condition(Span({'x': 2, 'y': [-1, 1]}), outlet), 'gamma_out': Condition(location=Span({'x': 2, 'y': [-1, 1]}), function=outlet),
'gamma_in': Condition(Span({'x': -2, 'y': [-1, 1]}), inlet), 'gamma_in': Condition(location=Span({'x': -2, 'y': [-1, 1]}), function=inlet),
'D': Condition(Span({'x': [-2, 2], 'y': [-1, 1]}), [momentum, continuity]), 'D1': Condition(location=Span({'x': [-2, 2], 'y': [-1, 1]}), function=momentum),
'D2': Condition(location=Span({'x': [-2, 2], 'y': [-1, 1]}), function=continuity),
} }
# conditions = {
# 'gamma_top': Condition(location=Span({'x': [-2, 2], 'y': 1}), function=wall),
# 'gamma_bot': Condition(location=Span({'x': [-2, 2], 'y': -1}), function=wall),
# 'gamma_out': Condition(location=Span({'x': 2, 'y': [-1, 1]}), function=outlet),
# 'gamma_in': Condition(location=Span({'x': -2, 'y': [-1, 1]}), function=inlet),
# 'D': Condition(location=Span({'x': [-2, 2], 'y': [-1, 1]}), function=[momentum, continuity]),
# }

View File

@@ -1,38 +1,10 @@
import argparse import argparse
import torch
from torch.nn import Softplus from torch.nn import Softplus
from pina.problem import SpatialProblem
from pina.operators import grad
from pina.model import FeedForward from pina.model import FeedForward
from pina import Condition, Span, Plotter, PINN from pina import Plotter, PINN
from problems.first_order_ode import FirstOrderODE
class FirstOrderODE(SpatialProblem):
x_rng = [0, 5]
output_variables = ['y']
spatial_domain = Span({'x': x_rng})
def ode(input_, output_):
y = output_
x = input_
return grad(y, x) + y - x
def fixed(input_, output_):
exp_value = 1.
return output_ - exp_value
def solution(self, input_):
x = input_
return x - 1.0 + 2*torch.exp(-x)
conditions = {
'bc': Condition(Span({'x': x_rng[0]}), fixed),
'dd': Condition(Span({'x': x_rng}), ode),
}
truth_solution = solution
if __name__ == "__main__": if __name__ == "__main__":
@@ -44,6 +16,7 @@ if __name__ == "__main__":
parser.add_argument("id_run", help="number of run", type=int) parser.add_argument("id_run", help="number of run", type=int)
args = parser.parse_args() args = parser.parse_args()
# define Problem + Model + PINN
problem = FirstOrderODE() problem = FirstOrderODE()
model = FeedForward( model = FeedForward(
layers=[4]*2, layers=[4]*2,
@@ -51,7 +24,6 @@ if __name__ == "__main__":
input_variables=problem.input_variables, input_variables=problem.input_variables,
func=Softplus, func=Softplus,
) )
pinn = PINN(problem, model, lr=0.03, error_norm='mse', regularizer=0) pinn = PINN(problem, model, lr=0.03, error_norm='mse', regularizer=0)
if args.s: if args.s:

View File

@@ -1,83 +1,84 @@
import argparse # import argparse
import numpy as np # import numpy as np
import torch # import torch
from torch.nn import Softplus # from torch.nn import Softplus
from pina import PINN, LabelTensor, Plotter # from pina import PINN, LabelTensor, Plotter
from pina.model import MultiFeedForward # from pina.model import MultiFeedForward
from problems.parametric_elliptic_optimal_control_alpha_variable import ( # from problems.parametric_elliptic_optimal_control_alpha_variable import (
ParametricEllipticOptimalControl) # ParametricEllipticOptimalControl)
class myFeature(torch.nn.Module): # class myFeature(torch.nn.Module):
""" # """
Feature: sin(x) # Feature: sin(x)
""" # """
def __init__(self): # def __init__(self):
super(myFeature, self).__init__() # super(myFeature, self).__init__()
def forward(self, x): # def forward(self, x):
t = (-x.extract(['x1'])**2+1) * (-x.extract(['x2'])**2+1) # t = (-x.extract(['x1'])**2+1) * (-x.extract(['x2'])**2+1)
return LabelTensor(t, ['k0']) # return LabelTensor(t, ['k0'])
class CustomMultiDFF(MultiFeedForward): # class CustomMultiDFF(MultiFeedForward):
def __init__(self, dff_dict): # def __init__(self, dff_dict):
super().__init__(dff_dict) # super().__init__(dff_dict)
def forward(self, x): # def forward(self, x):
out = self.uu(x) # out = self.uu(x)
p = LabelTensor((out.extract(['u_param']) * x.extract(['alpha'])), ['p']) # p = LabelTensor((out.extract(['u_param']) * x.extract(['alpha'])), ['p'])
return out.append(p) # return out.append(p)
if __name__ == "__main__": # if __name__ == "__main__":
parser = argparse.ArgumentParser(description="Run PINA") # parser = argparse.ArgumentParser(description="Run PINA")
group = parser.add_mutually_exclusive_group(required=True) # group = parser.add_mutually_exclusive_group(required=True)
group.add_argument("-s", "-save", action="store_true") # group.add_argument("-s", "-save", action="store_true")
group.add_argument("-l", "-load", action="store_true") # group.add_argument("-l", "-load", action="store_true")
args = parser.parse_args() # args = parser.parse_args()
opc = ParametricEllipticOptimalControl() # opc = ParametricEllipticOptimalControl()
model = CustomMultiDFF( # model = CustomMultiDFF(
{ # {
'uu': { # 'uu': {
'input_variables': ['x1', 'x2', 'mu', 'alpha'], # 'input_variables': ['x1', 'x2', 'mu', 'alpha'],
'output_variables': ['u_param', 'y'], # 'output_variables': ['u_param', 'y'],
'layers': [40, 40, 20], # 'layers': [40, 40, 20],
'func': Softplus, # 'func': Softplus,
'extra_features': [myFeature()], # 'extra_features': [myFeature()],
}, # },
} # }
) # )
pinn = PINN( # pinn = PINN(
opc, # opc,
model, # model,
lr=0.002, # lr=0.002,
error_norm='mse', # error_norm='mse',
regularizer=1e-8) # regularizer=1e-8)
if args.s: # if args.s:
pinn.span_pts( # pinn.span_pts(
{'variables': ['x1', 'x2'], 'mode': 'random', 'n': 100}, # {'variables': ['x1', 'x2'], 'mode': 'random', 'n': 100},
{'variables': ['mu', 'alpha'], 'mode': 'grid', 'n': 5}, # {'variables': ['mu', 'alpha'], 'mode': 'grid', 'n': 5},
locations=['D']) # locations=['D'])
pinn.span_pts( # pinn.span_pts(
{'variables': ['x1', 'x2'], 'mode': 'grid', 'n': 20}, # {'variables': ['x1', 'x2'], 'mode': 'grid', 'n': 20},
{'variables': ['mu', 'alpha'], 'mode': 'grid', 'n': 5}, # {'variables': ['mu', 'alpha'], 'mode': 'grid', 'n': 5},
locations=['gamma1', 'gamma2', 'gamma3', 'gamma4']) # locations=['gamma1', 'gamma2', 'gamma3', 'gamma4'])
pinn.train(1000, 20) # pinn.train(1000, 20)
pinn.save_state('pina.ocp') # pinn.save_state('pina.ocp')
else: # else:
pinn.load_state('pina.ocp') # pinn.load_state('pina.ocp')
plotter = Plotter() # plotter = Plotter()
plotter.plot(pinn, components='y', fixed_variables={'alpha': 0.01, 'mu': 1.0}) # plotter.plot(pinn, components='y', fixed_variables={'alpha': 0.01, 'mu': 1.0})
plotter.plot(pinn, components='u_param', fixed_variables={'alpha': 0.01, 'mu': 1.0}) # plotter.plot(pinn, components='u_param', fixed_variables={'alpha': 0.01, 'mu': 1.0})
plotter.plot(pinn, components='p', fixed_variables={'alpha': 0.01, 'mu': 1.0}) # plotter.plot(pinn, components='p', fixed_variables={'alpha': 0.01, 'mu': 1.0})
raise NotImplementedError('not available problem at the moment...')

View File

@@ -37,7 +37,9 @@ if __name__ == "__main__":
if args.s: if args.s:
pinn.span_pts(200, 'grid', locations=['gamma_top', 'gamma_bot', 'gamma_in', 'gamma_out']) pinn.span_pts(200, 'grid', locations=['gamma_top', 'gamma_bot', 'gamma_in', 'gamma_out'])
pinn.span_pts(2000, 'random', locations=['D']) # pinn.span_pts(2000, 'random', locations=['D'])
pinn.span_pts(2000, 'random', locations=['D1'])
pinn.span_pts(2000, 'random', locations=['D2'])
pinn.train(10000, 100) pinn.train(10000, 100)
with open('stokes_history_{}.txt'.format(args.id_run), 'w') as file_: with open('stokes_history_{}.txt'.format(args.id_run), 'w') as file_:
for i, losses in pinn.history_loss.items(): for i, losses in pinn.history_loss.items():