// GenTensors - A simple tensor library
// Copyright (C) 2012-2026 Eric Bechet
//
// See the LICENSE file for license information and contributions.
// Please report all bugs and problems to <bechet@cadxfem.org>.

// check gentensor ops
#include "genTensors.h"

int main(int argc,char**argv)
{
  double mat0=7.0;
  const double mat1[]={1.0,5.0,3.0};

  const double mat2[]=
  {
    1.0,0.0,-1.0,
    0.0,2.0,0.0,
    -1.0,0.0,3.0
  };

  const double mat3[]=
  {
    1.0,0.0,-1.0,
    0.0,2.0,0.0,
    -1.0,0.0,3.0,

    0.0,0.0,0.0,
    0.0,2.0,0.0,
    0.0,0.0,0.0,

    3.0,0.0,-10.0,
    0.0,4.0,0.0,
    -10.0,0.0,5.0
  };

  const double mat4[]=
  {
    1.0,0.0,-1.0,
    0.0,2.0,0.0,
    -1.0,0.0,3.0,

    0.0,0.0,0.0,
    0.0,2.0,0.0,
    0.0,0.0,0.0,

    3.0,0.0,-10.0,
    0.0,4.0,0.0,
    -10.0,0.0,5.0,


    1.0,0.0,-1.0,
    0.0,2.0,0.0,
    -1.0,0.0,3.0,

    0.0,0.0,0.0,
    0.0,2.0,0.0,
    0.0,0.0,0.0,

    3.0,0.0,-10.0,
    0.0,4.0,0.0,
    -10.0,0.0,5.0,


    1.0,0.0,-1.0,
    0.0,2.0,0.0,
    -1.0,0.0,3.0,

    0.0,0.0,0.0,
    0.0,2.0,0.0,
    0.0,0.0,0.0,

    3.0,0.0,-10.0,
    0.0,4.0,0.0,
    -10.0,0.0,5.0
  };


  genTensor0<double> t0(mat0);
  genTensor1<double> t1(mat1);
  genTensor2<double> t2(mat2);
  genTensor3<double> t3(mat3);
  genTensor4<double> t4(mat4);

  genTensor0<double> t00=t0*t0;
  genTensor0<double> t11=t1*t1;
  genTensor0<double> t22=t2*t2;
  genTensor0<double> t33=t3*t3;
  genTensor0<double> t44=t4*t4;
  
  genTensor1<double> t10=t1*t0;
  genTensor1<double> t01=t0*t1;
  genTensor1<double> t12=t1*t2;
  genTensor1<double> t21=t2*t1;
  genTensor1<double> t23=t2*t3;
  genTensor1<double> t32=t3*t2;
  genTensor1<double> t34=t3*t4;
  genTensor1<double> t43=t4*t3;

  genTensor2<double> t20=t2*t0;
  genTensor2<double> t02=t0*t2;
  genTensor2<double> t31=t3*t1;
  genTensor2<double> t13=t1*t3;
  genTensor2<double> t42=t4*t2;
  genTensor2<double> t24=t2*t4;

  genTensor3<double> t30=t3*t0;
  genTensor3<double> t03=t0*t3;
  genTensor3<double> t41=t4*t1;
  genTensor3<double> t14=t1*t4;

  genTensor4<double> t40=t4*t0;
  genTensor4<double> t04=t0*t4;

  std::cout << std::endl ;
  std::cout << std::endl ;
  std::cout.precision(17);
  std::cout << std::showpos;
  std::cout << std::scientific;

  std::cout <<  norm(t0) << " " <<  norm(t1) << " " <<  norm(t2) << " " <<  norm(t3) << " " <<  norm(t4) << std::endl;
  std::cout << norm(t00) << " " << norm(t11) << " " << norm(t22) << " " << norm(t33) << " " << norm(t44) << std::endl;
  std::cout << norm(t10) << " " << norm(t01) << " " << norm(t21) << " " << norm(t12) << " " << norm(t32) << " " << norm(t23) << " " << norm(t43) <<  " " << norm(t34) << std::endl;
  std::cout << norm(t20) << " " << norm(t02) << " " << norm(t31) << " " << norm(t13) << " " << norm(t42) << " " << norm(t24) << std::endl;
  std::cout << norm(t30) << " " << norm(t03) << " " << norm(t41) << " " << norm(t14) << std::endl;
  std::cout << norm(t40) << " " << norm(t04)  << std::endl;

  if (fabs(norm(t0) - 7.00000000000000000e+00)> 1e-12) return 1;
  if (fabs(norm(t1) - 5.91607978309961613e+00)> 1e-12) return 1;
  if (fabs(norm(t2) - 4.00000000000000000e+00)> 1e-12) return 1;
  if (fabs(norm(t3) - 1.64316767251549827e+01)> 1e-12) return 1;
  if (fabs(norm(t4) - 2.84604989415154144e+01)> 1e-12) return 1;

  if (fabs(norm(t00) - 49. ) >1e-12) return 1;
  if (fabs(norm(t11) - 35. ) >1e-12) return 1;
  if (fabs(norm(t22) - 16. ) >1e-12) return 1;
  if (fabs(norm(t33) - 65. ) >1e-12) return 1;
  if (fabs(norm(t44) - 1.  ) >1e-12) return 1;

  if (fabs(norm(t10) - 4.14125584816973102e+01)> 1e-12) return 1;
  if (fabs(norm(t01) - 4.14125584816973102e+01)> 1e-12) return 1;
  if (fabs(norm(t21) - 1.29614813968157208e+01)> 1e-12) return 1;
  if (fabs(norm(t12) - 1.29614813968157208e+01)> 1e-12) return 1;
  if (fabs(norm(t32) - 4.88671668914824622e+01)> 1e-12) return 1;
  if (fabs(norm(t23) - 3.76563407675254638e+01)> 1e-12) return 1;
  if (fabs(norm(t43) - 1.12583302491977022e+02)> 1e-12) return 1;
  if (fabs(norm(t34) - 9.13290753265355306e+01)> 1e-12) return 1;

  if (fabs(norm(t20) - 2.80000000000000000e+01)> 1e-12) return 1;
  if (fabs(norm(t02) - 2.80000000000000000e+01)> 1e-12) return 1;
  if (fabs(norm(t31) - 3.77094152699296075e+01)> 1e-12) return 1;
  if (fabs(norm(t13) - 5.40555270069583145e+01)> 1e-12) return 1;
  if (fabs(norm(t42) - 8.46404158779953093e+01)> 1e-12) return 1;
  if (fabs(norm(t24) - 3.28633534503099654e+01)> 1e-12) return 1;

  if (fabs(norm(t30) - 1.15021737076084889e+02)> 1e-12) return 1;
  if (fabs(norm(t03) - 1.15021737076084889e+02)> 1e-12) return 1;
  if (fabs(norm(t41) - 6.53146231712317302e+01)> 1e-12) return 1;
  if (fabs(norm(t14) - 1.47885090526394862e+02)> 1e-12) return 1;

  if (fabs(norm(t40) - 1.99223492590607890e+02)> 1e-12) return 1;
  if (fabs(norm(t04) - 1.99223492590607890e+02)> 1e-12) return 1;

  return 0; 
}