How to assign resolved shear stress for every slipsystem in Crystal plasticity

2

There are several issues that I can see by quickly skimming through it.

  • Edit tau looks good to me. You should just change the projection and it should work
  • Quadrature elements are finicky and should be used only when necessary. In any case, you will have a mismatch of quadrature points when trying to do a L2 projection. Your u is in lagrange space, wherereas tau isn’t.
  • Edit: tau is in Z which has nSlip.
    Try starting with a simpler example to get acquainted with FEniCS and then move your way up. This seems fairly involved at beginner level.

Some general code changes to make it more readable :slight_smile:


# Not advisable to set it globally unless you know it's accurate enough for your form
parameters["form_compiler"]["quadrature_degree"] = deg_stress 

def C_int(d):
    I = np.eye(d)
    II = np.einsum("ik,jl->ijkl", I, I) + np.einsum("il,jk->ijkl", I, I)
    JJ = np.einsum("ij,kl->ijkl", I, I)
    return as_tensor(mu*II + lamb*JJ)

# slip plane normals
sno_ = np.array(( 
        ((1, 1, 1),
        (1, 1, 1),
        (1, 1, 1),
        (-1, -1, 1),
        (-1, -1, 1),
        (-1, -1, 1),
        (-1, 1, 1),
        (-1, 1, 1),
        (-1, 1, 1),
        (1, -1, 1),
        (1,-1, 1),
        (1, -1, 1)) 
)) * 1./3**0.5

sdo_ = np.array((
             (0, 1, -1),
             (-1, 0, 1),
             (1, -1, 0),
             (0, -1, -1),
             (1, 0, 1),
             (-1, 1, 0),
             (0, 1, -1),
             (1, 0, 1),
             (-1, -1, 0),
             (0, -1, -1),
             (-1, 0, 1),
             (1, 1, 0)
    )) * 1.0/np.sqrt(2.0)

# get Resolved shear stress on all Slip systems 
def get_tau(J,S):
    for i in range(0,nSlip):
        mi, ni = sdo_[i], sno_[i]
        tau.sub(i).assign( project(1/J * dot(as_vector(mi), dot(S, as_vector(ni))) , Z0) ) # this should change
    return tau
3 Likes