-- FreeFem++ v4.14 (mer. 06 mars 2024 16:59:04 CET - git v4.14-1-g2b2052ae)
   file : lapack_test.edp
 Load: lg_fem lg_mesh lg_mesh3 eigenvalue 
    1 : //  lapack_test.edp
    2 : //
    3 : //  Licensing:
    4 : //
    5 : //    This code is distributed under the MIT license.
    6 : //
    7 : //  Modified:
    8 : //
    9 : //    12 November 2020
   10 : //
   11 : cout << "\n";
   12 : cout << "lapack_test()\n";
   13 : cout << "  FreeFem++ version.\n";
   14 : cout << "  Test the use of lapack functions.\n";
   15 : 
   16 : load "lapack" Add lapack interface ...
   17 : {
   18 : int n=5;
   19 : real[int,int] A(n,n);
   20 : real[int,int] A1(n,n);
   21 : real[int,int] B(n,n);
   22 : 
   23 : for(int i=0;i<n;++i)
   24 : {
   25 :   for(int j=0;j<n;++j)
   26 :   {
   27 :     A(i,j)= (i==j) ? n+1 : 1;
   28 :   }
   29 :  A(0,n-1)=-100;
   30 : }
   31 : cout << A << endl;
   32 : A1=A^-1;
   33 : cout << A1 << endl;
   34 : 
   35 : complex[int] vp(n);
   36 : complex[int,int] VP(n,n),KK(n,n);
   37 : 
   38 : int nn= dgeev(A,vp,VP);
   39 : cout << "  vp = " <<vp << endl;
   40 : cout << " VP = " << VP << endl;
   41 : KK =0.;
   42 : for(int i=0;i<n;++i)
   43 :   for(int j=0;j<n;++j)
   44 :     for(int k=0;k<n;++k)
   45 :       KK(i,j) += (A(i,k) - vp[j]* real(i==k) ) *VP(k,j);
   46 : cout <<" KK " <<  KK << endl; 
   47 : 
   48 : B=0;
   49 : for(int i=0;i<n;++i)
   50 :   for(int j=0;j<n;++j)
   51 :     for(int k=0;k<n;++k)
   52 :       B(i,j) +=A(i,k)*A1(k,j);
   53 : cout << B << endl;
   54 : B = A*A1; // version 3.13
   55 : cout << B << " ++++++  "<<  endl;
   56 : // A1+A^-1;  attention ne marche pas
   57 : 
   58 : inv(A1);
   59 : cout << A1 << endl;
   60 : }
   61 : {
   62 :   cout << " *****************  complex VP \n";
   63 : int n=5;
   64 : complex[int,int] A(n,n),A1(n,n),B(n,n);
   65 : for(int i=0;i<n;++i)
   66 : for(int j=0;j<n;++j)
   67 :   A(i,j)= (i==j) ? n+1 : 1;
   68 : A(0,n-1)=-100i;
   69 : cout << A << endl;
   70 : A1=A^-1;
   71 : cout << A1 << endl;
   72 : 
   73 :  complex[int] vp(n);
   74 :  complex[int,int] VP(n,n),KK(n,n);
   75 : 
   76 :  int nn= zgeev(A,vp,VP);
   77 :  cout << "  vp = " <<vp << endl;
   78 :  cout << " VP = " << VP << endl;
   79 :  KK =0.;
   80 :  for(int i=0;i<n;++i)
   81 :    for(int j=0;j<n;++j)
   82 :      for(int k=0;k<n;++k)
   83 :        KK(i,j) += (A(i,k) - vp[j]* real(i==k) ) *VP(k,j);
   84 :  cout <<" KK " <<  KK << endl; 
   85 :  assert( KK.linfty < 1e-10);
   86 : B=0;
   87 : for(int i=0;i<n;++i)
   88 :   for(int j=0;j<n;++j)
   89 :     for(int k=0;k<n;++k)
   90 :       B(i,j) +=A(i,k)*A1(k,j);
   91 : cout << B << endl;
   92 : // A1+A^-1;  attention ne marche pas
   93 : 
   94 : //inv(A1);
   95 : //cout << A1 << endl;
   96 :   cout << " *****************  fin --complex VP \n";
   97 : }
   98 : 
   99 : {
  100 :   int n=5;
  101 :   complex[int,int] A(n,n),A1(n,n),B(n,n);
  102 :   for(int i=0;i<n;++i)
  103 :     for(int j=0;j<n;++j)
  104 :       A(i,j)= (i==j) ? (n+1)*1i : 1+0i;
  105 :   cout << A << endl;
  106 :   A1=A^-1;
  107 :   cout << A1 << endl;
  108 : 
  109 :   B=0;
  110 :   for(int i=0;i<n;++i)
  111 :     for(int j=0;j<n;++j)
  112 :       for(int k=0;k<n;++k)
  113 : 	B(i,j) +=A(i,k)*A1(k,j);
  114 :   cout << B << endl;
  115 :   // A1+A^-1;  attention ne marche pas 
  116 : }
  117 : 
  118 : // FFCS - value for regression checks
  119 : real regtest=0;
  120 : 
  121 : {
  122 : 
  123 : int n=5;
  124 : real [int,int] A(n,n),  B(n,n),  C(n,n);
  125 : macro SETM(A,i,j,fij) { for(int i=0;i<A.n;++i) for(int j=0;j<A.m;++j) A(i,j)= (fij) ;} )  // 
  126 : SETM(C,i,j,1./(1+i+j))         { for(int i=0;i<C.n;++i) for(int j=0;j<C.m;++j) C(i,j)= (1./(1+i+j)) ;} ; 
  127 : SETM(B,i,j,i==j?2.:1./n)         { for(int i=0;i<B.n;++i) for(int j=0;j<B.m;++j) B(i,j)= (i==j?2.:1./n) ;} ;
  128 :  
  129 : A = B*C;
  130 : cout << A << " " << endl;
  131 : A = B + C;
  132 : real[int]  b(n),c(n);
  133 : real[int,int] AA=A;
  134 : AA=A;
  135 : real [int,int] At=A';
  136 : 
  137 : b = A*c;
  138 : 
  139 : At =A';
  140 : {
  141 : 
  142 : real [int,int] A=[[0,-1,2],[4,11,2],[0,-1,2],[4,11,2]];
  143 : real [int,int] B=[[3,-1],[1,2],[6,1]];
  144 : real [int,int] E=[[11,0],[35,20],[11,0],[35,20]];
  145 : real [int,int]  At=A';
  146 : real [int,int] C(A.n,B.m);
  147 : C=A*B;
  148 : cout << " C = " <<  C << endl;
  149 : cout << " E = " <<  E << endl;
  150 : 
  151 : C -= E;
  152 : assert( C.linfty < 1e-10);
  153 : C = A*B;
  154 : C -=E; 
  155 : assert( C.linfty < 1e-10);
  156 : }
  157 : {
  158 : 
  159 : complex [int,int] A=[[0,-1,2],[4,11,2],[0,-1,2],[4,11,2]];
  160 : complex [int,int] B=[[3,-1],[1,2],[6,1]];
  161 : complex [int,int] E=[[11,0],[35,20],[11,0],[35,20]];
  162 : 
  163 : complex [int,int] C(A.n,B.m);
  164 : complex [int,int]  At=A';
  165 : cout << "At = " << At << endl;
  166 : C=A*B;
  167 : cout << " C = " <<  C << endl;
  168 : cout << " E = " <<  E << endl;
  169 : 
  170 : C -= E;
  171 : assert( C.linfty < 1e-10);
  172 : C = A*B;
  173 : C -=E; 
  174 : assert( C.linfty < 1e-10);
  175 : 
  176 : // FFCS - value for regression checks
  177 : regtest=C.linfty;
  178 : }
  179 : }
  180 : //
  181 : //  Terminate.
  182 : //
  183 : cout << "\n";
  184 : cout << "lapack_test():\n";
  185 : cout << "  Normal end of execution.\n";
  186 : 
  187 : exit ( 0 );
  188 :  sizestack + 1024 =3696  ( 2672 )


lapack_test()
  FreeFem++ version.
  Test the use of lapack functions.
5 5	
	   6   1   1   1 -100
	   1   6   1   1   1
	   1   1   6   1   1
	   1   1   1   6   1
	   1   1   1   1   6
	
5 5	
	 0.05960264901 -0.140397351 -0.140397351 -0.140397351 1.063576159
	 -0.006622516556 0.1933774834 -0.006622516556 -0.006622516556 -0.140397351
	 -0.006622516556 -0.006622516556 0.1933774834 -0.006622516556 -0.140397351
	 -0.006622516556 -0.006622516556 -0.006622516556 0.1933774834 -0.140397351
	 -0.006622516556 -0.006622516556 -0.006622516556 -0.006622516556 0.05960264901
	
  vp = 5	
	(7.5,9.733961167)	(7.5,-9.733961167)	(5,0)	(5,0)	(5,0)
	
 VP = 5 5	
	 (0.9799979794,0) (0.9799979794,0) (0.8623065827,0) (0.8470625579,0) (0.6172702636,0)
	 (-0.01515460793,-0.09834291005) (-0.01515460793,0.09834291005) (-0.3623938326,0) (-0.1130263287,0) (-0.4503101658,0)
	 (-0.01515460793,-0.09834291005) (-0.01515460793,0.09834291005) (-0.2415394678,0) (-0.3546593284,0) (0.3649953345,0)
	 (-0.01515460793,-0.09834291005) (-0.01515460793,0.09834291005) (-0.2583732822,0) (-0.3793769009,0) (-0.5319554323,0)
	 (-0.01515460793,-0.09834291005) (-0.01515460793,0.09834291005) (-2.151389752e-17,0) (-3.612078844e-17,0) (-4.47886478e-17,0)
	
 KK 5 5	
	 (1.998401444e-15,0) (1.998401444e-15,0) (1.318722483e-15,0) (3.834123449e-15,0) (2.369441034e-15,0)
	 (-3.885780586e-16,3.053113318e-16) (-3.885780586e-16,-3.053113318e-16) (2.560418586e-16,0) (1.859238165e-16,0) (7.323674694e-16,0)
	 (-4.440892099e-16,1.942890293e-16) (-4.440892099e-16,-1.942890293e-16) (1.450195562e-16,0) (1.859238165e-16,0) (-1.48807858e-15,0)
	 (-3.885780586e-16,1.804112415e-16) (-3.885780586e-16,-1.804112415e-16) (1.450195562e-16,0) (1.859238165e-16,0) (9.544120744e-16,0)
	 (-4.440892099e-16,1.665334537e-16) (-4.440892099e-16,-1.665334537e-16) (-7.702504875e-17,0) (1.859238165e-16,0) (-4.888778577e-16,0)
	
5 5	
	   1   0   0 1.110223025e-16   0
	 -2.255140519e-17   1 -1.908195824e-17 -2.168404345e-17 1.665334537e-16
	 -3.295974604e-17 -6.245004514e-17   1 -2.168404345e-17 3.60822483e-16
	 -3.469446952e-17 -6.418476861e-17 -6.938893904e-18   1 1.387778781e-16
	 -1.387778781e-17 -5.551115123e-17 1.387778781e-17 6.938893904e-18   1
	
5 5	
	   1 3.469446952e-17 -2.775557562e-17 5.898059818e-17 3.330669074e-16
	 -2.255140519e-17   1 -1.908195824e-17 -2.168404345e-17 1.665334537e-16
	 -3.122502257e-17 -6.071532166e-17   1 -2.168404345e-17 3.053113318e-16
	 -3.816391647e-17 -6.765421556e-17 -1.040834086e-17   1 1.387778781e-16
	 -1.387778781e-17 -5.204170428e-17 1.387778781e-17 8.67361738e-18   1
	 ++++++  
5 5	
	   6   1   1   1 -100
	   1   6   1   1   1
	   1   1   6   1   1
	   1   1   1   6   1
	   1   1   1   1   6
	
 *****************  complex VP 
5 5	
	 (6,0) (1,0) (1,0) (1,0) (0,-100)
	 (1,0) (6,0) (1,0) (1,0) (1,0)
	 (1,0) (1,0) (6,0) (1,0) (1,0)
	 (1,0) (1,0) (1,0) (6,0) (1,0)
	 (1,0) (1,0) (1,0) (1,0) (6,0)
	
5 5	
	 (0.0364256805,-0.07142290294) (-0.1635743195,-0.07142290294) (-0.1635743195,-0.07142290294) (-0.1635743195,-0.07142290294) (1.272168875,0.6428061265)
	 (-0.004047297834,0.007935878105) (0.1959527022,0.007935878105) (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (-0.1635743195,-0.07142290294)
	 (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (0.1959527022,0.007935878105) (-0.004047297834,0.007935878105) (-0.1635743195,-0.07142290294)
	 (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (0.1959527022,0.007935878105) (-0.1635743195,-0.07142290294)
	 (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (-0.004047297834,0.007935878105) (0.0364256805,-0.07142290294)
	
  vp = 5	
	(0.2409464229,6.887950264)	(14.75905358,-6.887950264)	(5,0)	(5,0)	(5,0)
	
 VP = 5 5	
	 (0.9842714065,0) (0.9760761115,0) (-0.14023168,4.062520713e-17) (0.8290735144,0) (0.8232205804,0)
	 (-0.06943417571,-0.05460169237) (0.06973372763,0.08340301773) (0.8387536883,0) (-0.5019051387,-8.128818563e-17) (-0.2002154408,-1.850251993e-17)
	 (-0.06943417571,-0.05460169237) (0.06973372763,0.08340301773) (-0.477442181,-3.224086188e-18) (-0.2236207034,-6.492363346e-17) (-0.1015573969,-2.938635518e-17)
	 (-0.06943417571,-0.05460169237) (0.06973372763,0.08340301773) (-0.2210798273,-2.241172349e-18) (-0.1035476723,-5.301017147e-17) (-0.5214477427,-2.938635518e-17)
	 (-0.06943417571,-0.05460169237) (0.06973372763,0.08340301773) (-7.632281524e-17,-1.393534548e-18) (4.112970128e-17,-4.273633023e-17) (3.48282728e-17,-2.938635518e-17)
	
 KK 5 5	
	 (-2.664535259e-15,0) (-5.329070518e-15,8.881784197e-16) (1.937134526e-16,7.667441472e-15) (-5.161811443e-15,-4.312192118e-15) (-3.271702425e-15,-3.56010251e-15)
	 (-3.330669074e-16,-2.775557562e-16) (-4.440892099e-16,8.326672685e-17) (-2.407791167e-15,3.376641405e-17) (3.741966087e-16,-2.419583208e-16) (-2.982386346e-16,-1.066615855e-16)
	 (-3.330669074e-16,-2.775557562e-16) (0,8.326672685e-17) (1.089411361e-15,3.376641405e-17) (1.243964281e-16,-2.419583208e-16) (-2.982386346e-16,-1.066615855e-16)
	 (-3.330669074e-16,-2.775557562e-16) (0,6.938893904e-17) (3.955219702e-16,3.376641405e-17) (2.725191347e-17,-2.419583208e-16) (-2.982386346e-16,-1.066615855e-16)
	 (-3.330669074e-16,-2.775557562e-16) (0,5.551115123e-17) (-1.873451177e-16,3.376641405e-17) (-6.989260119e-17,-2.419583208e-16) (-2.982386346e-16,-1.066615855e-16)
	
5 5	
	 (1,0) (1.110223025e-16,0) (0,0) (0,5.551115123e-17) (-8.881784197e-16,-4.440892099e-16)
	 (-4.250072516e-17,-8.67361738e-18) (1,-1.040834086e-17) (-1.474514955e-17,-8.67361738e-18) (-3.729655473e-17,8.67361738e-18) (-2.498001805e-16,-6.938893904e-17)
	 (-3.556183126e-17,-2.255140519e-17) (-3.035766083e-17,-1.734723476e-17) (1,-2.255140519e-17) (-9.540979118e-18,-1.734723476e-18) (0,-9.714451465e-17)
	 (-2.949029909e-17,-5.204170428e-18) (-2.428612866e-17,-3.469446952e-18) (-1.561251128e-17,-5.204170428e-18) (1,-5.204170428e-18) (2.775557562e-17,-1.249000903e-16)
	 (-4.510281038e-17,-1.387778781e-17) (-4.857225733e-17,-2.081668171e-17) (-3.122502257e-17,-1.387778781e-17) (-2.775557562e-17,-1.387778781e-17) (1,-1.110223025e-16)
	
 *****************  fin --complex VP 
5 5	
	 (0,6) (1,0) (1,0) (1,0) (1,0)
	 (1,0) (0,6) (1,0) (1,0) (1,0)
	 (1,0) (1,0) (0,6) (1,0) (1,0)
	 (1,0) (1,0) (1,0) (0,6) (1,0)
	 (1,0) (1,0) (1,0) (1,0) (0,6)
	
5 5	
	 (-0.006237006237,-0.1528066528) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356)
	 (0.02079002079,0.009355509356) (-0.006237006237,-0.1528066528) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356)
	 (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (-0.006237006237,-0.1528066528) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356)
	 (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (-0.006237006237,-0.1528066528) (0.02079002079,0.009355509356)
	 (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (0.02079002079,0.009355509356) (-0.006237006237,-0.1528066528)
	
5 5	
	 (1,-1.734723476e-18) (3.469446952e-18,-6.938893904e-18) (-3.469446952e-18,-1.734723476e-18) (6.938893904e-18,-1.734723476e-18) (6.071532166e-18,-2.775557562e-17)
	 (6.938893904e-18,-2.949029909e-17) (1,6.938893904e-18) (1.040834086e-17,-1.734723476e-18) (1.387778781e-17,-5.724587471e-17) (6.071532166e-18,-2.775557562e-17)
	 (0,-1.561251128e-17) (3.469446952e-18,6.938893904e-18) (1,1.214306433e-17) (6.938893904e-18,-1.734723476e-18) (-8.67361738e-19,-2.775557562e-17)
	 (-1.387778781e-17,2.602085214e-17) (3.469446952e-18,5.204170428e-17) (-3.469446952e-18,-1.734723476e-18) (1,5.204170428e-18) (-4.33680869e-18,0)
	 (0,2.775557562e-17) (6.938893904e-18,2.775557562e-17) (1.387778781e-17,-2.775557562e-17) (6.938893904e-18,0) (1,6.938893904e-18)
	
5 5	
	 2.256666667 1.19 0.8185714286 0.6269047619 0.5091269841
	 1.356666667 0.89 0.6685714286 0.5369047619 0.4491269841
	 1.056666667 0.74 0.5785714286 0.4769047619 0.4062698413
	 0.9066666667 0.65 0.5185714286 0.434047619 0.3741269841
	 0.8166666667 0.59 0.4757142857 0.4019047619 0.3491269841
	 
 C = 4 2	
	  11   0
	  35  20
	  11   0
	  35  20
	
 E = 4 2	
	  11   0
	  35  20
	  11   0
	  35  20
	
At = 3 4	
	 (0,0) (4,0) (0,0) (4,0)
	 (-1,0) (11,0) (-1,0) (11,0)
	 (2,0) (2,0) (2,0) (2,0)
	
 C = 4 2	
	 (11,0) (0,0)
	 (35,0) (20,0)
	 (11,0) (0,0)
	 (35,0) (20,0)
	
 E = 4 2	
	 (11,0) (0,0)
	 (35,0) (20,0)
	 (11,0) (0,0)
	 (35,0) (20,0)
	

lapack_test():
  Normal end of execution.
  current line = 187
exit(0)
 err code 0 ,  mpirank 0
 CodeAlloc : nb ptr  4825,  size :548448 mpirank: 0
Ok: Normal End