#include "sparselizardbase.h"

using namespace mathop;

void sparselizard(void)
{
    formulation elasticity;
    int bottom = 1, top = 2, bulk = 3;
    double lambda, mu;
    
    setaxisymmetry();
    mesh mymesh("axi-lizard.msh");

    field u("h1xyz");
    u.setorder(bulk, 2);
    u.setconstraint(bottom);
  
    parameter E, nu;
    E|bulk = 100e3; nu|bulk = 0.3;
  
    elasticity += integral(bulk, predefinedelasticity(dof(u), tf(u), E, nu));
    elasticity += integral(top, array1x3(0,-100,0)*tf(u));

    elasticity.generate();
    solve(elasticity);

    u.write(bulk, "u.vtk");
}

int main(void)
{	
    SlepcInitialize(0,{},0,0);

    sparselizard();

    SlepcFinalize();

    return 0;
}

