


GETVISCOSITY2D - returns viscosity from a tensor, using Glen's flow law
From a string tensor and a material object, return viscosity, using Glen's flow law.
2*B
mu= -------------------------------------------------------------------
2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n]
where mu is the viscotiy, B the flow law parameter , (u,v) the velocity
vector, and n the flow law exponent.
If epsilon is NULL, it means this is the first time Emg is being run, and we
return 10^14, initial viscosity.
This routine is appropriate for 2d elements.
Usage:
mu=GetViscosity2d(matice,epsilon)

0001 function mu=GetViscosity2d(matice,epsilon) 0002 %GETVISCOSITY2D - returns viscosity from a tensor, using Glen's flow law 0003 % 0004 % From a string tensor and a material object, return viscosity, using Glen's flow law. 0005 % 2*B 0006 % mu= ------------------------------------------------------------------- 0007 % 2[ exx^2+eyy^2+exx*eyy+exy^2+exz^2+eyz^2 ]^[(n-1)/2n] 0008 % 0009 % where mu is the viscotiy, B the flow law parameter , (u,v) the velocity 0010 % vector, and n the flow law exponent. 0011 % 0012 % If epsilon is NULL, it means this is the first time Emg is being run, and we 0013 % return 10^14, initial viscosity. 0014 % 0015 % This routine is appropriate for 2d elements. 0016 % 0017 % Usage: 0018 % mu=GetViscosity2d(matice,epsilon) 0019 0020 eps0=10^-27; 0021 if matice.n==1, 0022 %Viscous behaviour. mu=*B; 0023 mu=matice.B; 0024 else 0025 %Non linear creep 0026 if norm(epsilon,2)~=0, 0027 0028 B=matice.B; 0029 n=matice.n; 0030 0031 exx=epsilon(1); 0032 eyy=epsilon(2); 0033 exy=epsilon(3); 0034 0035 %Build viscosity: mu=2*B/(2*A^e): 0036 A=exx^2+eyy^2+exy^2+exx*eyy+eps0^2; 0037 e=(n-1)/2/n; 0038 mu=2*B/(2*A^e); 0039 0040 else 0041 mu=10^14; 0042 end 0043 end