1 | function vector_average=DepthAverage(md,vector)
|
---|
2 | %DEPTHAVERAGE - computes depth average of 3d vector using the trapezoidal rule, and returns the value on 2d mesh.
|
---|
3 | %
|
---|
4 | % Usage:
|
---|
5 | % vector_average=DepthAverage(md,vector);
|
---|
6 | %
|
---|
7 | % Example:
|
---|
8 | % vel_bar=DepthAverage(md,md.initialization.vel);
|
---|
9 |
|
---|
10 | %check that the model given in input is 3d
|
---|
11 | if ~strcmp(md.mesh.elementtype(),'Penta');
|
---|
12 | error('DepthAverage error message: the model given in input must be 3d')
|
---|
13 | end
|
---|
14 |
|
---|
15 | %nods data
|
---|
16 | if (length(vector)==md.mesh.numberofvertices),
|
---|
17 | vector_average=zeros(md.mesh.numberofvertices2d,1);
|
---|
18 | for i=1:md.mesh.numberoflayers-1,
|
---|
19 | vector_average=vector_average+(project2d(md,vector,i)+project2d(md,vector,i+1))/2.*(project2d(md,md.mesh.z,i+1)-project2d(md,md.mesh.z,i));
|
---|
20 | end
|
---|
21 | vector_average=vector_average./project2d(md,md.geometry.thickness,1);
|
---|
22 |
|
---|
23 | %element data
|
---|
24 | elseif (length(vector)==md.mesh.numberofelements),
|
---|
25 | vector_average=zeros(md.mesh.numberofelements2d,1);
|
---|
26 | for i=1:md.mesh.numberoflayers-1,
|
---|
27 | vertices_dz = (project2d(md,md.mesh.z,i+1)-project2d(md,md.mesh.z,i));
|
---|
28 | elements_dz = mean(vertices_dz(md.mesh.elements2d),2);
|
---|
29 | vector_average = vector_average+project2d(md,vector,i).*elements_dz;
|
---|
30 | %vector_average=vector_average+project2d(md,vector,i).*(project2d(md,md.mesh.z,i+1)-project2d(md,md.mesh.z,i));
|
---|
31 | end
|
---|
32 | vertices_thickness = project2d(md,md.geometry.thickness,1);
|
---|
33 | elements_thickness = mean(vertices_thickness(md.mesh.elements2d),2);
|
---|
34 | vector_average = vector_average./elements_thickness;
|
---|
35 | %vector_average=vector_average./project2d(md,md.geometry.thickness,1);
|
---|
36 |
|
---|
37 | else
|
---|
38 | error('vector size not supported yet');
|
---|
39 | end
|
---|