femwind/vertical_test.m

   1 %test for vertical_sweeps
   2
   3 format compact
   4 % test for vertical_sweeps
   5 nel=[5,4,3]
   6 %nel=[1 1 1]
   7 n=nel+1
   8 h = [1,1,1]
   9 expand=1.3
  10 A = diag([1 1 1])
  11 lambda=[]
  12 params=[]
  13 u0=[];
  14
  15 X = regular_mesh(nel,h,expand);
  16 % X = add_terrain_to_mesh(X, 'hill', 'shift', 0.1)
  17 X = add_terrain_to_mesh(X, 'hill', 'squash', 0.1)  % more thorough testing
  18
  19 % plot_mesh_3d(X)
  20
  21 % assemble sparse system matrix
  22 [K,F,~] = sparse_assembly(A,X,u0,lambda,params);
  23
  24 nn = size(F,1);
  25 %xr=rand(n1,n2,n3);
  26 x = rand(nn,1);
  27
  28 K_1=nd_assembly(A,X,lambda,params);
  29 [n1,n2,n3,m1,m2,m3]=size(K_1);
  30 K27 = reshape(K_1,n1,n2,n3,m1*m2*m3);  % make to n x 27 nd format
  31 K14 = ndt_convert(K27,14);
  32
  33 % test same results for ndt_mult from matlab and fortran
  34 if exist('fortran/sweeps_test.exe')
  35     disp('testing if same result in fortran')
  36     err=sweeps_fortran(K,K14,n1,n2,n3,F,X,x)
  37     if abs(err)<1e-6
  38         fprintf('error %g OK\n',err)
  39     else
  40         error(sprintf('error %g too large',err))
  41     end
  42 else
  43     warning('fortran/sweeps_test.exe not available')
  44 end