I find that the matlab official example "deflectio​n-of-a-pie​zoelectric​-actuator" cannot be correct

While teaching piezoelectrical materials in my course at the Politecnico di BARI, I have explored the Matlab example "deflection-of-a-piezoelectric-actuator" which can be found in the official matlab web site at https://it.mathworks.com/help/pde/ug/deflection-of-a-piezoelectric-actuator.html

The result of the example itself may be correct because it is checked against some reference solution.

However this is when the voltage across the bilayer actuator is simply linear across the y-direction because uniform voltage is imposed on the top surface and on the bottom one. However, when I asked the students to reproduce other conditions, like for example voltage only across part of the top surface and on the bottom one, or across the lateral surfaces, the resulting electrical field is very weird and very unlikely to be correct.

After all, although the problem is strictly speaking coupled between the electrical field and the mechanical one, I would not expect the electrical field to be very different from the classical electrical or thermal problem in the same domain, with imposed voltage (or temperature). And this problem, solved with Matlab or with other means, shows much smoother results.

The treatment of matrices is so involute that it is hard to check if the formulation is correct.

But when we checked the solution using ANSYS, even the coupled piezoelectrical problem shows smooth results, which differ from this example.

So be aware that there must be some error in the formulation!

Prof. Michele Ciavarella

Politecnico di BARI

[e-mail]

Answer 1

[Ravi Kumar]

Hello Prof. Ciavarella,

Can you provide a bit more details, reproduction code and expected results, of the modified setup to investigate this issue?

Regards,

Ravi

[michele ciavarella]

thanks Ravi

the short comment is that if you change the B.C. applying voltage NOT on the top and bottom surfaces (1, 2 respectively), but on lateral ones (6,7 with 100 volts, and 3,4 with 0 volts), you don't find any reasonable result --- the proper solution should have voltage varying smoothly from 100 to 0 volts moving from left to rigth). Just change

Set the voltage (solution component 3) on the top edge to V.

voltTop = applyBoundaryCondition(model,'mixed','Edge',1,...

'u',V,...

'EquationIndex',3);

Set the voltage (solution component 3) on the bottom edge to zero.

voltBot = applyBoundaryCondition(model,'mixed','Edge',2,...

'u',0,...

'EquationIndex',3);

into

voltTop1 = applyBoundaryCondition(model,'mixed','Edge',6,...

'u',V,...

'EquationIndex',3);

voltTop2 = applyBoundaryCondition(model,'mixed','Edge',7,...

'u',V,...

'EquationIndex',3);

voltBot1 = applyBoundaryCondition(model,'mixed','Edge',3,...

'u',0,...

'EquationIndex',3);

voltBot2 = applyBoundaryCondition(model,'mixed','Edge',4,...

'u',0,...

'EquationIndex',3);

and see the results.

I am asking my students to provide you an example in ANSYS with the solution which is indeed qualitatively different.

mike

[Ravi Kumar]

Hello Prof. Ciavarella,

Thanks for the detailed reproduction steps. If I follow your steps exactly, i.e., replace original voltTop and voltBot with four new BCs and don't do any other change, then I do see voltage is zero across and not vary linearly from left to right. This is expected as left the clampLeft BC which basically overwrites the new voltTop1 and voltTop2. Hence, if I comment the clampLeft BC then I get the expected linear variation of the voltage as shown in the figure below.

Let me know if I have missed anything.

Regards,

Ravi

[michele ciavarella]

That's interesting, although I don't see this should be correct. There seems to be interference between the mechanical b.c. and the electrical ones. The ANSYS solution with the clamped end shows linear variation in voltage too. Also, I am not sure this is the only problem. The students created other b.c., like voltage across only part of the top and bottom surfaces, and the voltage also appears weird --- like it doesn't transmits well in the x-direction.

[Ravi Kumar]

Hello Prof. Ciavarella,

If you don't want to overwrite voltage BCs with structural BCs, then you need to specify all three equations dirichlet BC in once call as:

voltTop1 = applyBoundaryCondition(model,'dirichlet','Edge',6,...
                                               'u',[0,0,V],...
                                               'EquationIndex',1:3);
voltTop2 = applyBoundaryCondition(model,'dirichlet','Edge',7,...
                                               'u',[0,0,V],...
                                               'EquationIndex',1:3);                                          
                                           

I would be happy to look at the issue you are facing with voltage applied only on a portion of top and bottom surfaces if you can provide the reproduction steps.

Regards,

Ravi

[michele ciavarella]

dear Ravi

my students have further looked at this problem, and they have concluded that they do not know how to apply b.c. of both mechanical and electrical type on the same boundary. If they apply mechanical b.c. on one edge (Dirichlet or mixed, it doesn't matter), and electrical on another, with two different command lines, the code seems to work ok.

Hence, the point is how to apply for example clamping plus a given voltage, on the same edge?

regards, mike

[Ravi Kumar]

Hello Prof. Ciavarella,

My reply on 16 May 2019 at 22:46 precisely shows what you want to do. There are three equations in the setup, first two are structural equation and third is for potential. Clamped structural BC corresponds to setting first two dependent variables to zero. Hence, the following two BCs sets zero displacement and a voltage on edge 6 and 7.

voltLeft1 = applyBoundaryCondition(model,'dirichlet','Edge',6,...
                                               'u',[0,0,V],...
                                               'EquationIndex',1:3);
voltLeft2 = applyBoundaryCondition(model,'dirichlet','Edge',7,...
                                               'u',[0,0,V],...
                                               'EquationIndex',1:3);  
                                           

You can also write the above two statements as one statement with both 6 and 7 in the Edge option as:

voltLeft1n2 = applyBoundaryCondition(model,'dirichlet','Edge',[6,7],...
                                               'u',[0,0,V],...
                                               'EquationIndex',1:3);
                                           

Regards,

Ravi

[michele ciavarella]

Thank you very much. This issue of how to introduce BCs explains the problem.