SDT-piezo         Contents     Functions         Previous Next     PDF Index

3.5  Example 1 : Static response of a piezoelectric patch

3.5.1  Static response of a patch in extension mode

In this very simple example, the electric field and the strains are all constant, so that the electric potential and the displacement field are linear. The example is treated analytically in section section 2.4. It is therefore possible to obtain an exact solution using a single volumic 8-node finite element (with linear shape functions, the nodal unknowns being the displacements in x,y and z and the electric potential φ). Consider a piezoelectric patch whose dimensions and material properties are given in Table 3.1. The material properties correspond to the material SONOX_P502_iso in m_piezo.


PropertyValue
b10 mm
w10 mm
h2 mm
E54 GPa
ν0.41
d31=d32-185  10−12 pC/N (or m/V)
d33440  10−12 pC/N (or m/V)
d15=d24560  10−12 pC/N (or m/V)
є33T22T11T1850 є0
є08.854 10−12 Fm−1
Table 3.1: Geometrical and material properties of the piezoelectric patch

We first produce the mesh, associate the material properties and define the electrodes with
d_piezo('TutoPatch-s1') . The default material is SONOX_P502_iso. The number of elements in the x, y and z directions are given by nx,ny and nz.

% See full example as MATLAB code in d_piezo('ScriptTutoPatch')
t_avc('DefineStyles');

%% Step 1 Build mesh - Define electrodes
% Meshing script can be viewed with sdtweb d_piezo('MeshPatch')
model=d_piezo('MeshPatch lx=1e-2 ly=1e-2 h=2e-3 nx=1 ny=1 nz=1');
% Define electrodes
model=p_piezo('ElectrodeMPC Top -ground',model,'z==2e-3');
model=p_piezo('ElectrodeMPC Bottom -Input "Free patch"',model,'z==0');

The information about the nodes associated to each electrode can be obtained through the following call:

p_piezo('TabInfo',model)

The material can be changed for example to PIC_255 with the following call, and the full set of mechanical, piezoelectric and permittivity matrices can be obtained in order to check consistency with the datasheet (d_piezo('TutoPatch-s2') ):

%% Step 2 Define material properties
model.pl=m_piezo('dbval 1 -elas 2 PIC_255');
p_piezo('TabDD',model) % creates the table with full set of matrices

The next step consists in defining the boundary conditions and load case using
d_piezo('TutoPatch-s3') . We consider here two cases, the first one where the patch is free to expand, and the second one where it is mechanically constrained (all mechanical degrees of freedom are equal to 0).

%% Step 3 Compute static response
% to avoid rigid body mode
model=stack_set(model,'info','Freq',10);
def=fe_simul('dfrf',model); def.lab={'Free patch, axial'};
def.fun=[0 1]; def=feutil('rmfield',def,'data','LabFcn');

% Append mechanically constrained structure
%    can't call fe_simul because no free DOF
%    see code with sdtweb d_piezo('scriptFullConstrain')
def=d_piezo('scriptFullConstrain',model,def);
def.lab{2}='Constrained patch, axial';

We can look at the deformed shape, and plot the electric field for both cases.
(d_piezo('TutoPatch-s4')

%% Step 4 Visualize deformed shape
cf=feplot(model,def);
% Electric field representation
p_piezo('viewElec EltSel "matid1" DefLen 20e-4 reset',cf);
fecom('colormap',[1 0 0]);fecom('undef line');iimouse('resetview');
cf.mdl.name='E-field'; % Sets figure name
t_avc('SetStyle',cf); feplot(cf);

Figure 3.4: Vizualisation of the electric field and deformed shape for the free patch under unit voltage excitation

For the free patch deformed shape, we compute the mean strains from which d31, d32 and d33 are deduced. The values are found to be equal to the analytical values used in the model. Note that the parameters of the constitutive equations can be recovered using (d_piezo('TutoPatch-s5') ):

%% Step 5 : check constitutive law
% Decompose constitutive law
CC=p_piezo('viewdd -struct',cf); %

where the fields of CC are self-explanatory. The parameters which are not directly defined are computed from the equations presented in Section section 2.1.

% Display and compute mean strains
a=p_piezo('viewstrain -curve -mean',cf); % Strain S
fprintf('Relation between mean strain on free structure and d_3i\n');
E3=a.Y(9,1); disp({'E3 mean' a.Y(9,1)  1/2e-3 'E3 analytic'})

disp([a.X{1}(1:3) num2cell([a.Y(1:3,1)/E3 CC.d(3,1:3)']) ...
     {'d_31';'d_32';'d_33'}])

For the constrained patch, we compute the mean stress from which we can compute the e31, e32 and e33 values which are found to be equal to the analytical values used in the model:

% Display and compute mean stresses
b=p_piezo('viewstress -curve -mean',cf); % Stress T
fprintf('Relation between mean stress on pure electric and e_3i\n');
disp([b.X{1}(1:3) num2cell([b.Y(1:3,2)/-E3 CC.e(3,1:3)']) ...
     {'e_31';'e_32';'e_33'}])

% Mean stress/strain
disp([b.X{1} num2cell(b.Y(:,2)) num2cell(a.Y(:,1)) a.X{1}])

We can also compute the charge and the charge density (in pC/m2) accumulated on the electrodes, and compare with the analytical values (d_piezo('TutoPatch-s6') ):

%% Step 6 Check capacitance values
% Theoretical values of Capacitance and charge density - free patch
CT=CC.epst_r(3,3)*8.854e-12*1e-2*1e-2/2e-3; %% Capacitance - free patch
CdensT=CC.epst_r(3,3)*8.854e-12/2e-3*1e12; %% charge density - free patch

% Theoretical values of Capacitance and charge density - constrained patch
CS=CC.epss_r(3,3)*8.854e-12*1e-2*1e-2/2e-3; %% Capacitance - free patch
CdensS=CC.epss_r(3,3)*8.854e-12/2e-3*1e12; %% charge density - free patch

% Represent charge density (C/S) value on the electrodes
%   - compare with analytical values
cut=p_piezo('electrodeviewcharge',cf,struct('EltSel','matid 1'));
b=fe_caseg('stressobserve',cut,cf.def);b=reshape(b.Y,[],2);
disp([{'','CdensT','CdensS'};{'Numeric';'Theoretical'} ...
  num2cell([mean(abs(b));CdensT CdensS])])
  iimouse('zoom reset');

% Compute the value of the total charge (from reaction at electrical dof)
% Compare with analytical values
p_piezo('electrodeTotal',cf) %
disp('Theoretical values of capacitance')
disp([{'CT';'CS'} num2cell([CT;CS])])
cf.mdl.name='Charge'; % Sets figure name

t_avc('SetStyle',cf); feplot(cf);
fecom('ch 1')


   
Figure 3.5: Vizualisation of the total charge on the electrodes for the unconstrained and constrained patch under unit voltage excitation

The results clearly show the very large difference of charge density between the two cases (free patch or constrained patch).
For this simple static example, a finer mesh can be used, but it does not lead to more accurate results (this can be done by changing the values in the call to d_piezo('mesh') for example:

% Build mesh with refinement
model=d_piezo('MeshPatch lx=1e-2 ly=1e-2 h=2e-3 nx=5 ny=5 nz=2');
% Now a model with quadratic elements
model=d_piezo('MeshPatch lx=1e-2 ly=1e-2 h=2e-3 Quad');

3.5.2  Static response of a patch in shear mode

As for the patch in extension, as the fields are also uniform (see section 2.5), the problem can be modelled with a single 8-node element. The patch is meshed and then the poling is aligned with the −y axis by performing a rotation of 90o around the x-axis (d_piezo('TutoPatchShear-s1') ).

t_avc('DefineStyles');

% See full example as MATLAB code in d_piezo('ScriptTutoShearPatch')
%% Step 1 Build mesh and define electrodes
%Meshing script can be viewed with sdtweb d_piezo('MeshPatch')
model=d_piezo('MeshPatch lx=1e-2 ly=1e-2 h=2e-3 nx=1 ny=1 nz=1');

% Define electrodes
model=p_piezo('ElectrodeMPC Top -ground',model,'z==2e-3');
model=p_piezo('ElectrodeMPC Bottom -Input "Free patch"',model,'z==0');

% Rotate basis to align poling direction with y (-90 deg around x)
model.bas=basis('rotate',[],'rx=-90',1); %create local basis with id=1
model=feutil('setpro 1 COORDM=1',model); % assign basis with id=1 to pro=1

Then the response is computed both for the free case and the fully constrained case. The deformed shape for the free case is shown in Figure 3.6 together with the applied electric field:
(d_piezo('TutoShearPatch-s2') )


Figure 3.6: Deformed shape of a piezoelectric patch poled in the plane with an electric field applied in the out-of-plane direction

%% Step 2 Compute static response
% to avoid rigid body mode
model=stack_set(model,'info','Freq',10);
def=fe_simul('dfrf',model); def.lab={'Free patch, shear'};
def.fun=[0 1]; def=feutil('rmfield',def,'data','LabFcn');

% Append mechanically constrained structure
%    can't call fe_simul because no free DOF
%    see code with sdtweb d_piezo('scriptFullConstrain')
def=d_piezo('scriptFullConstrain',model,def);
def.lab{2}='Constrained patch, shear';

(d_piezo('TutoShearPatch-s3') )

%% Step 3 Vizualise deformed shape
cf=feplot(model,def); fecom('undef line');

% Electric field representation
p_piezo('viewElec EltSel "matid1" DefLen 20e-4 reset',cf);
cf.mdl.name='E-field'; % Sets figure name
t_avc('SetStyle',cf); feplot(cf);

The mean of shear strain and stress is evaluated and compared to the d24 piezo coefficient. Note that the mean values are computed in the global yz axis for which a negative strain corresponds to a positive strain in the local 23 axis. Finally the capacitance is evaluated and compared to the theoretical values, showing a perfect agreement, and demonstrating the difference with the extension case for CS.
(d_piezo('TutoShearPatch-s4') )

%% Step 4 : Check constitutive law
% Decompose constitutive law
CC=p_piezo('viewdd -struct',cf); %

% Display and compute mean strains
a=p_piezo('viewstrain -curve -mean',cf); % Strain S
fprintf('Relation between mean strain on free structure and d_24\n');
E3=a.Y(9,1); disp({'E3 mean' a.Y(9,1)  1/2e-3 'E3 analytic'})

disp([a.X{1}(4) num2cell([a.Y(4,1)/E3 CC.d(2,4)']) ...
     {'d_24'}])

% Display and compute mean stresses
b=p_piezo('viewstress -curve -mean',cf); % Stress T
fprintf('Relation between mean stress on pure electric and e_24 \n');
disp([b.X{1}(4) num2cell([b.Y(4,2)/-E3 CC.e(2,4)']) ...
     {'e_24'}])
% Mean stress/strain
disp([b.X{1} num2cell(b.Y(:,2)) num2cell(a.Y(:,1)) a.X{1}])

% Theoretical values of Capacitance and charge density - free patch
CT=CC.epst_r(2,2)*8.854e-12*1e-2*1e-2/2e-3; %% Capacitance - free patch
CdensT=CC.epst_r(2,2)*8.854e-12/2e-3*1e12; %% charge density - free patch

% Theoretical values of Capacitance and charge density - constrained patch
CS=CC.epss_r(2,2)*8.854e-12*1e-2*1e-2/2e-3; %% Capacitance - constrained patch
CdensS=CC.epss_r(2,2)*8.854e-12/2e-3*1e12; %% charge density - constrained patch

% Represent charge density (C/S) value on the electrodes
%   - compare with analytical values
cut=p_piezo('electrodeviewcharge',cf,struct('EltSel','matid 1'));
b=fe_caseg('stressobserve',cut,cf.def);b=reshape(b.Y,[],2);
disp([{'','CdensT','CdensS'};{'Numeric';'Theoretical'} ...
  num2cell([mean(abs(b));CdensT CdensS])])
    iimouse('zoom reset');

(d_piezo('TutoShearPatch-s5') )

%% Step 5 Check capacitance
% Compute the value of the total charge (from reaction at electrical dof)
% Ccompare with analytical values
p_piezo('electrodeTotal',cf) %
disp('Theoretical values of capacitance')
disp([{'CT';'CS'} num2cell([CT;CS])])

©1991-2024 by SDTools
Previous Up Next