Ask Question

Name:
Title:
Your Question:

Answer Question

Name:
Your Answer:
User Submitted Source Code!


Description:
  file.pas
Language: PASCAL
Code:
unit USparseRealLinearSolverLapacke;



interface

uses
  Classes, SysUtils, USparseLinearSolver, UTypes;

type TSparseRealLinearSolverLapacke = class (TSparseRealLinearSolver)
  public
    procedure Solve;override;
end;

implementation

uses ulapacke,ctypes;

procedure TSparseRealLinearSolverLapacke.solve;
var
  info:lapack_int;
  ipiv:array of cint;
  a:array of cdouble; //Matrix will be there
  b: TDoubleArray; //Vector for RHS and Solution
  n,ii:cint;
begin
  checkSizes;
  n:=prMatrix.MatrixOrder;
  SetLength(a,n*n);
  SetLength(ipiv,n);
  for ii:=0 to n*n-1 do
    a[ii]:=0.0;
  for ii:=0 to prMatrix.getLength()-1 do
    with prMatrix.getElementFromIndex(ii) do
      a[i*n+j]:=a[i*n+j]+value;

  SetLength(b, Length(prRHSVector.value));
  Move(prRHSVector.Value[0], b[0], n * SizeOf(prRHSVector.Value[0]));

  info := LAPACKE_dsysv(
    LAPACK_COL_MAJOR,cchar('U'),n,1,a,n,ipiv,b,n);

  LapackAssert('dsysv', info);

  prSolVector.value := b;
end;

end.
          
Comments: