# pdenonlin

Solve nonlinear elliptic PDE problem

## Syntax

• ``u = pdenonlin(model,c,a,f)``
example
• ``u = pdenonlin(b,p,e,t,c,a,f)``
example
• ``u = pdenonlin(___,Name,Value)``
example
• ``````[u,res] = pdenonlin(___)``````

## Description

example

````u = pdenonlin(model,c,a,f)` solves the nonlinear PDE$-\nabla \cdot \left(c\nabla u\right)+au=f,$with geometry, boundary conditions, and finite element mesh in `model`, and coefficients `c`, `a`, and `f`. In this context, "nonlinear" means some coefficient in `c`, `a`, or `f` depends on the solution `u` or its gradient. If the PDE is a system of equations (`model.PDESystemSize` > 1), then `pdenonlin` solves the system of equations$-\nabla \cdot \left(c\otimes \nabla u\right)+au=f.$```

example

````u = pdenonlin(b,p,e,t,c,a,f)` solves the PDE with boundary conditions `b`, and finite element mesh (`p`,`e`,`t`).```

example

````u = pdenonlin(___,Name,Value)`, for any previous arguments, modifies the solution process with `Name`, `Value` pairs.```
``````[u,res] = pdenonlin(___)``` also returns the norm of the Newton step residuals `res`.```

## Examples

collapse all

### Minimal Surface Problem

Solve a minimal surface problem. Because this problem has a nonlinear c coefficient, use `pdenonlin` to solve it.

Create a model and include circular geometry using the built-in `circleg` function.

```model = createpde; geometryFromEdges(model,@circleg);```

Set the coefficients.

```a = 0; f = 0; c = '1./sqrt(1+ux.^2+uy.^2)';```

Set a Dirichlet boundary condition with value x2.

```boundaryfun = @(region,state)region.x.^2; applyBoundaryCondition(model,'Edge',1:model.Geometry.NumEdges,... 'u',boundaryfun,'Vectorized','on');```

Generate a mesh and solve the problem.

```generateMesh(model,'Hmax',0.1); u = pdenonlin(model,c,a,f); pdeplot(model,'xydata',u,'zdata',u);```

### Minimal Surface Problem Using Legacy Approach

Solve the minimal surface problem using the legacy approach for creating boundary conditions and geometry.

Create the geometry using the built-in `circleg` function. Plot the geometry to see the edge labels.

```g = @circleg; pdegplot(g,'EdgeLabels','on') axis equal```

Create Dirichlet boundary conditions with value x2. Create the following file and save it on your MATLAB® path. For details of this approach, see Boundary Conditions by Writing Functions.

```function [qmatrix,gmatrix,hmatrix,rmatrix] = pdex2bound(p,e,u,time) ne = size(e,2); % number of edges qmatrix = zeros(1,ne); gmatrix = qmatrix; hmatrix = zeros(1,2*ne); rmatrix = hmatrix; for k = 1:ne x1 = p(1,e(1,k)); % x at first point in segment x2 = p(1,e(2,k)); % x at second point in segment xm = (x1 + x2)/2; % x at segment midpoint y1 = p(2,e(1,k)); % y at first point in segment y2 = p(2,e(2,k)); % y at second point in segment ym = (y1 + y2)/2; % y at segment midpoint switch e(5,k) case {1,2,3,4} hmatrix(k) = 1; hmatrix(k+ne) = 1; rmatrix(k) = x1^2; rmatrix(k+ne) = x2^2; end end```

Set the coefficients and boundary conditions.

```a = 0; f = 0; c = '1./sqrt(1+ux.^2+uy.^2)'; b = @pdex2bound;```

Generate a mesh and solve the problem.

```[p,e,t] = initmesh(g,'Hmax',0.1); u = pdenonlin(b,p,e,t,c,a,f); pdeplot(p,e,t,'xydata',u,'zdata',u);```

### Nonlinear Problem with 3-D Geometry

Solve a nonlinear 3-D problem with nontrivial geometry.

Import the geometry from the `BracketWithHole.stl` file. Plot the geometry and face labels.

```model = createpde(); importGeometry(model,'BracketWithHole.stl'); h = pdegplot(model,'FaceLabels','on'); h(1).FaceAlpha = 0.5;```

Set a Dirichlet boundary condition with value 1000 on the back face, which is face `3`. Set the large faces `1` and `5`, and also the circular face `9`, to have Neumann boundary conditions with value `g` = –10. Do not set boundary conditions on the other faces. Those faces default to Neumann boundary conditions with value `g` = 0.

```applyBoundaryCondition(model,'Face',3,'u',1000); applyBoundaryCondition(model,'Face',[1,5,9],'g',-10);```

Set the `c` coefficient to 1, `f` to 0.1, and `a` to the nonlinear value ```'0.1 + 0.001*u.^2'```.

```c = 1; f = 0.1; a = '0.1 + 0.001*u.^2';```

Generate the mesh and solve the PDE. Start from the initial guess `u0` = 1000, which matches the value you set on face `3`. Turn on the `Report` option to observe the convergence during the solution.

```generateMesh(model); u = pdenonlin(model,c,a,f,'U0',1000,'Report','on');```
```Iteration Residual Step size Jacobian: full 0 3.5985e-01 1 1.0182e-01 1.0000000 2 3.0318e-02 1.0000000 3 8.7499e-03 1.0000000 4 1.9023e-03 1.0000000 5 1.5581e-04 1.0000000 6 1.2424e-06 1.0000000```

Plot the solution on the geometry boundary.

`pdeplot3D(model,'colormapdata',u);`

## Input Arguments

collapse all

### `model` — PDE model`PDEModel` object

PDE model, specified as a `PDEModel` object.

Example: `model = createpde(1)`

### `c` — PDE coefficientscalar or matrix | character array | coefficient function

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. `c` represents the c coefficient in the scalar PDE

$-\nabla \cdot \left(c\nabla u\right)+au=f,$

or in the system of PDEs

$-\nabla \cdot \left(c\otimes \nabla u\right)+au=f.$

You can specify`c` in various ways, detailed in c Coefficient for Systems. See also Scalar PDE Coefficients, Specify Scalar PDE Coefficients in String Form, Specify 2-D Scalar Coefficients in Function Form, Specify 3-D PDE Coefficients in Function Form, and Coefficients for Systems of PDEs.

Example: `'cosh(x+y.^2)'`

Data Types: `double` | `char` | `function_handle`
Complex Number Support: Yes

### `a` — PDE coefficientscalar or matrix | character array | coefficient function

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. `a` represents the a coefficient in the scalar PDE

$-\nabla \cdot \left(c\nabla u\right)+au=f,$

or in the system of PDEs

$-\nabla \cdot \left(c\otimes \nabla u\right)+au=f.$

You can specify`a` in various ways, detailed in a or d Coefficient for Systems. See also Scalar PDE Coefficients, Specify Scalar PDE Coefficients in String Form, Specify 2-D Scalar Coefficients in Function Form, Specify 3-D PDE Coefficients in Function Form, and Coefficients for Systems of PDEs.

Example: `2*eye(3)`

Data Types: `double` | `char` | `function_handle`
Complex Number Support: Yes

### `f` — PDE coefficientscalar or matrix | character array | coefficient function

PDE coefficient, specified as a scalar or matrix, as a character array, or as a coefficient function. `f` represents the f coefficient in the scalar PDE

$-\nabla \cdot \left(c\nabla u\right)+au=f,$

or in the system of PDEs

$-\nabla \cdot \left(c\otimes \nabla u\right)+au=f.$

You can specify`f` in various ways, detailed in f Coefficient for Systems. See also Scalar PDE Coefficients, Specify Scalar PDE Coefficients in String Form, Specify 2-D Scalar Coefficients in Function Form, Specify 3-D PDE Coefficients in Function Form, and Coefficients for Systems of PDEs.

Example: `char('sin(x)';'cos(y)';'tan(z)')`

Data Types: `double` | `char` | `function_handle`
Complex Number Support: Yes

### `b` — Boundary conditionsboundary matrix | boundary file

Boundary conditions, specified as a boundary matrix or boundary file. Pass a boundary file as a function handle or as a string naming the file.

For more information on boundary conditions, see Forms of Boundary Condition Specification.

Example: `b = 'circleb1'` or equivalently ```b = @circleb1```

Data Types: `double` | `char` | `function_handle`

### `p` — Mesh nodesoutput of `initmesh` | output of `meshToPet`

Mesh nodes, specified as the output of `initmesh` or `meshToPet`. For the structure of a `p` matrix, see Mesh Data for [p,e,t] Triples: 2-D and Mesh Data for [p,e,t] Triples: 3-D.

Example: `[p,e,t] = initmesh(g)`

Data Types: `double`

### `e` — Mesh edgesoutput of `initmesh` | output of `meshToPet`

Mesh edges, specified as the output of `initmesh` or `meshToPet`. For the structure of `e`, see Mesh Data for [p,e,t] Triples: 2-D and Mesh Data for [p,e,t] Triples: 3-D.

Example: `[p,e,t] = initmesh(g)`

Data Types: `double`

### `t` — Mesh elementsoutput of `initmesh` | output of `meshToPet`

Mesh elements, specified as the output of `initmesh` or `meshToPet`. Mesh elements are the triangles or tetrahedra that form the finite element mesh. For the structure of a `t` matrix, see Mesh Data for [p,e,t] Triples: 2-D and Mesh Data for [p,e,t] Triples: 3-D.

Example: `[p,e,t] = initmesh(g)`

Data Types: `double`

### Name-Value Pair Arguments

Specify optional comma-separated pairs of `Name,Value` arguments. `Name` is the argument name and `Value` is the corresponding value. `Name` must appear inside single quotes (`' '`). You can specify several name and value pair arguments in any order as `Name1,Value1,...,NameN,ValueN`.

Example:

### `'Jacobian'` — Approximation of Jacobian`'full'` (3-D default) | `'fixed'` (2-D default) | `'lumped'`

Approximation of Jacobian, specified as `'full'`, `'fixed'`, or `'lumped'`.

• `'full'` means numerical evaluation of the full Jacobian based on the sparse version of the `numjac` function. 3-D geometry uses only `'full'`, any other specification yields an error.

• `'fixed'` specifies a fixed-point iteration matrix where the Jacobian is approximated by the stiffness matrix. This is the 2-D geometry default.

• `'lumped'` specifies a "lumped" approximation as described in Nonlinear Equations. This approximation is based on the numerical differentiation of the coefficients.

Example: `u = pdenonlin(model,c,a,f,'Jacobian','full')`

Data Types: `char`

### `'U0'` — Initial solution guess0 (default) | scalar | string | vector

Initial solution guess, specified as a scalar, string, or vector. For details, see Initial Conditions.

• A scalar specifies a constant initial condition for either a scalar or PDE system.

• For scalar problems, write a string using the same syntax as Specify Scalar PDE Coefficients in String Form.

• For systems of N equations, write a character array with N rows, where each row has the syntax of Specify Scalar PDE Coefficients in String Form.

• For systems of N equations, and a mesh with `Np` nodes, give a column vector with N*`Np` components. The nodes are either `model.Mesh.Nodes`, or the `p` data from `initmesh` or `meshToPet`. See Mesh Data.

The first `Np` elements contain the values of component 1, where the value of element `k` corresponds to node `p(k)`. The next `Np` points contain the values of component 2, etc. It can be convenient to first represent the initial conditions `u0` as an `Np`-by-`N` matrix, where the first column contains entries for component 1, the second column contains entries for component 2, etc. The final representation of the initial conditions is `u0(:)`.

Example: `u = pdenonlin(model,c,a,f,'U0','x.^2-y.^2')`

Data Types: `double` | `char`
Complex Number Support: Yes

### `'Tol'` — Residual size at termination1e-4 (default) | positive scalar

Residual size at termination, specified as a positive scalar. `pdenonlin` iterates until the residual size is less than `'Tol'`.

Example: `u = pdenonlin(model,c,a,f,'Tol',1e-6)`

Data Types: `double`

### `'MaxIter'` — Maximum number of Gauss-Newton iterations`25` (default) | positive integer

Maximum number of Gauss-Newton iterations, specified as a positive integer.

Example: `u = pdenonlin(model,c,a,f,'MaxIter',12)`

Data Types: `double`

### `'MinStep'` — Minimum damping of search direction`1/2^16` (default) | positive scalar

Minimum damping of search direction, specified as a positive scalar.

Example: `u = pdenonlin(model,c,a,f,'MinStep',1e-3)`

Data Types: `double`

### `'Report'` — Print convergence information`'off'` (default) | `'on'`

Print convergence information, specified as `'off'` or `'on'`.

Example: `u = pdenonlin(model,c,a,f,'Report','on')`

Data Types: `char`

### `'Norm'` — Residual norm`Inf` (default) | p value for Lp norm | `'energy'`

Residual norm, specified as the `p` value for Lp norm, or as the string `'energy'`. `p` can be any positive real value, `Inf`, or `-Inf`. The `p` norm of a vector `v` is `sum(abs(v)^p)^(1/p)`. See `norm`.

Example: `u = pdenonlin(model,c,a,f,'Norm',2)`

Data Types: `double` | `char`

## Output Arguments

collapse all

### `u` — PDE solutionvector

PDE solution, returned as a vector.

• If the PDE is scalar, meaning only one equation, then `u` is a column vector representing the solution u at each node in the mesh. `u(i)` is the solution at the `i`th column of `model.Mesh.Nodes` or the `i`th column of `p`.

• If the PDE is a system of N > 1 equations, then `u` is a column vector with N*`Np` elements, where `Np` is the number of nodes in the mesh. The first `Np` elements of `u` represent the solution of equation 1, then next `Np` elements represent the solution of equation 2, etc.

To obtain the solution at an arbitrary point in the geometry, use `pdeInterpolant`.

To plot the solution, use `pdeplot` for 2-D geometry, or see Plot 3-D Solutions.

### `res` — Norm of Newton step residualsscalar

Norm of Newton step residuals, returned as a scalar. For information about the algorithm, see Nonlinear Equations.

expand all

### Tips

• If the Newton iteration does not converge, `pdenonlin` displays the error message `Too many iterations` or ```Stepsize too small```.

• If the initial guess produces matrices containing `NaN` or `Inf` elements, `pdenonlin` displays the error message ```Unsuitable initial guess U0 (default: U0 = 0)```.

### Algorithms

The `pdenonlin` algorithm is a Gauss-Newton iteration scheme applied to the finite element matrices. For details, see Nonlinear Equations.