femwind/femwind_solve.m

   1 function [W,rate]=femwind_solve(A,X,U0,params)
   2     % assemble sparse system matrix
   3     nel=size(X{1})-1;
   4     lambda = zeros(prod(nel+1),1); % placeholder solution
   5     % F = f_assembly_fortran(A,X,U0,lambda,params);
   6     [K,F,~] = sparse_assembly(A,X,U0,lambda,params);
   7     if isfield(params,'femwind_fortran_test') && params.femwind_fortran_test
   8         F_f=read_array('F_f');
   9         err_f = big(F(:)-F_f(:))
  10     end
  11
  12
  13     % dirichlet boundary conditions
  14     [K,~]=apply_boundary_conditions(K,[],X);   % to K - once
  15     [~,F]=apply_boundary_conditions([],F,X);   % to F - every time
  16
  17     % solve the equations
  18     % [lambda,it] = sparse_solve(K,F,X,'s');
  19     [lambda,it,rate,XC] = sparse_solve(K,F,X,params);
  20     format long
  21     rate
  22
  23     % assemble final wind
  24     [~,~,W] = sparse_assembly(A,X,U0,lambda,params);
  25 end