Ax = 0 with 2 values of x known

I have a matrix in the form Ax = 0, and for x I know the first and last position. How can I solve this system?

Answer 1

[Riccardo Scorretti]

Basically, you are asking how to impose Dirichlet boundary conditions, right? If so, there are several methods:

• Substitution method: replace the first and last equation by: x = known_value. For instance

A(1,:) = 0 ; A(1,1) = 1; b(1) = known_value_1;

A(n,:) = 0 ; A(n,n) = 1; b(n) = known_value_2;

• Penalty method: "overwrite" the existing equation with the boundary condition, weighted with a huge value, for instance:

TGV = 1E30; % = Terribly Giant Value

A(1,1) = TGV ; b(1) = TGV*known_value_1;

A(n,n) = TGV ; b(n) = TGV*known_value_2;

Assuming that b is the known vector, then you can solve the linear system with an appropriate direct or iterative linear solver.

The penalty method is approximate, but it has the advantage that it is quite simple and (most importantly), it preserves eventual symmetry properties of the matrix. This is the way (I don't know if it is the only way) Dirichlet boundary conditions are implemented in FreeFem++ (see https://doc.freefem.org/references/types.html):

See also this discussion:

https://scicomp.stackexchange.com/questions/5072/how-to-properly-apply-non-homogeneous-dirichlet-boundary-conditions-with-fem

By the way, I recently posted a complete program where Dirichlet boundary conditions are imposed by using the first method. The program solves a PDE by using the Finite Difference method, but the way of imposing Dirichlet boundary conditions is exactly the same:

https://fr.mathworks.com/matlabcentral/answers/1683989-second-order-boundary-value-problem-numerical-solution?s_tid=srchtitle

[Torsten]

@Riccardo Scorretti

I wonder what the size of A in the original question is.

Your answer applies if it is (n-2) x n, isn't it ?

[Riccardo Scorretti]

No, A must (and usually is) a square matrix. In FEM assembly algorithm, usually Dirichlet boundary conditions are imposed by modifying the linear system, which normally has as much unknown as equation by construction.

The point is that some equations are replaced by Dirichlet boundary conditions: you must have a nxn linear system before and after imposing Dirichlet BCs. It sounds logical: equations impose constraints between the unknowns, hence if you want to impose Dirichlet boundary condtions you have to remove some of the already existing constraints, otherwise you will obtain an over-determined system.

I don't think things are much different with other methods (= Finite Differences, Finite Volumes, or whatever). By the way, it is what I did in my code, which implements a FD method.

[Torsten]

Will substituting the first and the last equation by

A(1,:) = 0 ; A(1,1) = 1; b(1) = known_value_1;
A(n,:) = 0 ; A(n,n) = 1; b(n) = known_value_2;

lead to a solution for which also the original first and last equation is satisfied ?

I see that this is not necessary if the matrix A stems from the discretization of a differential equation and you just insert the conditional equation for x(1) and x(n).

But was the set-up of the matrix A as (n x n) not wrong from the beginning, i.e. assuming the equations for x(1) and x(n) as if they were inner grid points ? Shouldn't it have been (n-2) x n (discretization scheme in inner grid points) + Dirichlet conditions in boundary points right from the beginning ?

[Riccardo Scorretti]

Well, you are absolutely right. Indeed the question didn't mention explicitly that the linear systems comes from a PDE: I extrapolated (perhaps I'm wrong) from the way the question is posed. Otherwise the problem is over/under-derermined, unless A is (n-2)xn. In this case it sounds logical to solve it in the sense of constrained least squares.

[kimfet]

Hello,

This worked for me! It was Dirichlett BC for a FEM problem with a Galerkin approach, to solve the Steady-State Convection-Diffusion equation in 1D. A pretty easy CFD problem, but working on the solution of the system was a bit tricky. It worked for me (I used the substitution method because it was easier for me to implement).

Thank's for your time!

Answer 2

[Torsten]

Use lsqlin for the problem

min : ||A*x||_2
s.c.
x(1) = known_value_1
x(n) = known_value_2

if A has n columns.

Answer 3

[Matt J]

b=-A(:,[1,end])*[xInitial;xFinal];
x=A(:,2:end-1)\b;

x=[xInitial, x', xFinal]';

[Torsten]

b = -(xInitial*A(:,1) + xFinal*A(:,end));
x = A(:,2:end-1)) \ b;
x = [xInitial ; x ; xFinal];

[Matt J]

Right.

Answer 4

[Fabio Freschi]

The post is old but the question arises many times, so I post here my solution. It is a generalization of @Matt J solution.

Suppose you have a matrix with unknowns partitioned in free x and dirichlet . And by chance the partition is such that the free variables come first

The free unknowns can be obtained by solving the first block row:

The second block row is of no interest, since is known.

In the general case when are not at the end of the unknown vector, but they are identified by the index vector idFix and the corresponding Dirichlet values xFix, you can write the function

function x = solveproblem(A,b,idFix,xFix)
% get number of unknowns
n = length(b);
% move to rhs the Dirichlet values
b = b-A(:,idfix)*xfix;
% get indices of free variables 
iFree = setdiff(1:n,idFix);
% reduce stiffness and rhs
A = A(:,iFree);
A = A(iFree,:);
b = b(iFree,:);

% solution
xRed = A\b;

% re-assemble the unknown vector
x(idFree) = xRed;
x(idFix) = xFix;
end

This method preserves the original properties of the original matrix (in particular if it's SPD). There is also an interesting way to include in this approach that a selection of the vector of unknowns has the same value, for example a floating potential of an equipotential object, but it goes beyond the original question