applications/solvers/solidMechanics/elasticOrthoAcpSolidFoam/calculateTraction.H

   1 {
   2     surfaceVectorField n = mesh.Sf()/mesh.magSf();
   3
   4 //     traction = (n&fvc::interpolate(sigma));
   5
   6     // surfaceTensorField sGradU =
   7     //     ((I - n*n)&fvc::interpolate(gradU));
   8
   9     // traction =
  10     //     (2*mu + lambda)*snGradU
  11     //   - (mu + lambda)*(snGradU&(I - n*n))
  12     //   + mu*(sGradU&n)
  13     //   + lambda*tr(sGradU&(I - n*n))*n;
  14
  15 //     traction =
  16 //       (2*mu + lambda)*fvc::snGrad(U)
  17 //     - (mu + lambda)*(n&sGradU)
  18 //     + mu*(sGradU&n)
  19 //     + lambda*tr(sGradU)*n;
  20
  21     // philipc
  22     // I am having trouble with back-calculation of interface tractions from solid interface
  23     // procedure (in multiMaterial.C), the tractions have quite large differences from each
  24     // side. Interpolating sigma is OK for now
  25     // traction = rheology.law().interfaceTraction(n, U, gradU, rheology.mu(), rheology.lambda());
  26 #   include "calculateEpsilonSigma.H"
  27     traction = (n&fvc::interpolate(sigma));
  28
  29     // forAll(traction.boundaryField(), patchi)
  30     //   {
  31     //  if (mesh.boundary()[patchi].type() == "cohesive")
  32     //    {
  33     //      forAll(traction.boundaryField()[patchi], facei)
  34     //        {
  35     //          Pout << "face " << facei << " with traction magnitude "
  36     //               << mag(traction.boundaryField()[patchi][facei])/1e6 << " MPa and traction "
  37     //               << traction.boundaryField()[patchi][facei]/1e6 << " MPa" << endl;
  38     //        }
  39     //    }
  40     //   }
  41 }