Domain dependent parameters in elasticity

General Discussion
1

Hi All,

I want to implement elasticity in domain [0,1]*[0,1] with known analytical solution, and I want to use different values of lambda and nu, in different halves of the domain and I tried it as follows:

real nu1 = 10;  
real lambda1 = 2*10^4; Lame constant

real nu2 = 20;  
real lambda2 = 10^4; 

func nu = (y>=0.5 )*(nu1 + 0*x + 0*y) + (y<0.5 )*(nu2+ 0*x + 0*y);
func lambda = (y>=0.5 )*(lambda1 + 0*x + 0*y) + (y<0.5 )*(lambda2+ 0*x + 0*y);

On doing this, I noticed that at y=0.5, the L2 norm of computed solution u and exact solution uex are very different and I do not get the convergence of my solution. When I use the same parameters nu and lambda, I do not have any problem with the elasticity solver (I get desired convergence rates).

Please let me know if I needed to change anything in defining nu and lambda differently in different parts of the domain.

Thank you!

2

In the documentation section 2.6 (thermal conduction problem), different diffusion (k) has been used in different parts of the domain. So, I thought it should work the same way!

1)Defining nu and lambda as functions of x and y.
2) Use it in defining right-side functions and in variational formulation.

I plotted my solution in the case of elasticity in the whole domain [0,1]*[0,1] and it is way off as compared to exact solution. While using same lambda and nu over entire domain, my solution matches with the exact solution.

//RHS for Elasticity 
func f1s = -nus*(2*u1struexx+u1strueyy+u2struexy)-lambdas*(u1struexx+u2strueyx); 
func f2s = -nus*(2*u2strueyy+u2struexx+u1strueyx)-lambdas*(u1struexy+u2strueyy); 

fespace Uh(Gh,P2); // element space for displacement
Uh u1s,u2s, v1s, v2s;

problem structure([u1s,u2s], [v1s,v2s]) =
		int2d(Gh) ( 
						2*nus*(
								dx(u1s)*dx(v1s) 
								+ 0.5*(dx(u2s) + dy(u1s))*(dx(v2s) + dy(v1s)) 
								+ dy(u2s)*dy(v2s)
								) // 2nu(D(u),D(v))
					+ lambdas*(dx(u1s) + dy(u2s))*(dx(v1s) + dy(v2s)) // lambda(div(u),div(v))	
					)
		-	int2d(Gh)(
						f1s*v1s + f2s*v2s 
					)

		+ on(2,3,4,5,u1s=u1strue,u2s=u2strue); // BC
3

@fb77 Sir, I would really appriciate if you could please add any thought on it? Thank you!