applications/solvers/solidMechanics/elasticIncrAcpSolidFoam/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 "calculateDEpsilonDSigma.H"
  27     traction = (n&fvc::interpolate(sigma+DSigma));
  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 }