How can I find Volume integration of velocity

Rakib Hossain (view profile)

on 18 May 2019
Latest activity Commented on by Walter Roberson

Walter Roberson (view profile)

on 19 May 2019
I am trying to perform volume integration of velocity. These data are random but within a certain range.
x is from -10 to 30, y and z is from -1 to 1 and U is from 0 to 1.
My code is as follows-
U= data(:,1); % (:,1) means all first columns
x= data(:,2); % (:,2) means all second columns
y = data(:,3); % (:,3) means all third columns
z = data(:,4); % (:,4) means all fourth columns
result=trapz(x,trapz(y,trapz(z,U,3),2),1); %% My question is here.
These yields an error showing, Length(x) must equal the length of the DIM'th dimendion of Y.
How can I cope with this?

R2011a

Walter Roberson (view profile)

on 18 May 2019

U = data(:,1); % (:,1) means all first columns
So U is a column vector
z = data(:,4); % (:,4) means all fourth columns
so z is a column vector with the same size as U
trapz(z,U,3)
requests that the values in z have a trapazoid summation done along the third dimension, using the content of U as the coordinates along the third dimension. However we know that z is a column vector, so its length along the third dimension is 1, which does not match the length of U.
If it happens that x, y, z form a regular cuboid with value U at each point... then you do not have enough U values.
So we have to deduce that you do not have a regular grid: you have scattered points. You cannot use trapz to meaningfully do volume integration in such a case.
Perhaps you could reshape your x, y, z to form a regular cuboid. If so then you would need to know the order they are in, and you would need to know the marginal lengths of the edges. You would then reshape, such as
marglen = [17 11 16];
X = reshape(x, marglen);
Y = reshape(y, marglen);
Z = reshape(z, marglen);
UU = reshape(U, marglen);
and once that was done it would become meaningful to start using trapz over the various dimensions,
trapz( squeeze(X(:,1,1)), trapz(squeeze(Y(1,:,1)), trapz(squeeze(Z(1,1,:)), UU, 3), 2), 1)
Here the intent of the squeeze() section is to extract the marginal unique values for each of the directions.

Show 1 older comment
Walter Roberson

Walter Roberson (view profile)

on 18 May 2019
In the special case where your x, y, and z are arranged as a regular cuboid, there are several different ways that the data can be arranged:
1. x changes most rapidly, then y within x, then z most slowly
2. x changes most rapidly, then z within x, then y most slowly
3. y changes most rapidly, then x within y, then z most slowly
4. y changes most rapidly, then z within y, then x most slowly
5. z changes most rapidly, then x within z, then y most slowly
6. z changes most rapidly, then y within z, then x most slowly
You have to examine the data to see what it looks like. For example does it look like
10 1.1 7.5
15 1.1 7.5
20 1.1 7.5
10 1.2 7.5
15 1.2 7.5
20 1.2 7.5
10 1.1 8.4
15 1.1 8.4
20 1.1 8.4
10 1.2 8.4
15 1.2 8.4
20 1.2 8.4
10 1.1 9.3
15 1.1 9.3
20 1.1 9.3
10 1.2 9.3
15 1.2 9.3
20 1.2 9.3
This is a case where the x varies most quickly, then y within x, then z most slowly.
The above example has 3 different x values, 2 different y values, and 3 different z values. That would give marglen of [3 2 3] -- number of different x, number of different y, number of different z.
The code would have to be altered a bit if the order were different, such as if the third column were the one varying most quickly, such as
10 1.1 7.5
10 1.1 8.4
10 1.1 9.3
10 1.2 7.5
10 1.2 8.4
10 1.2 9.3
15 1.1 7.5
and so on
Neither of these orders is more "right" then the other, so you would have to examine the way that your data is and change the details accordingly.
However, this only applies if the data can be arranged in a cuboid. If you just have a list of x, y, z points that are not regular, then you need to rethink your approach.
Rakib Hossain

Rakib Hossain (view profile)

on 19 May 2019
My updated code is now
U = data(:,1); % (:,1) means all first columns
x = data(:,2); % (:,2) means all second columns
y = data(:,3); % (:,3) means all tthird columns
z = data(:,4); % (:,4) means all fourth columns
x_sort= sort(x);
y_sort= sort(y);
z_sort= sort(z);
%distinct element counting:
countx=0;
county=0;
countz=0;
for i=1:length(x)-1
if(x_sort(i)~=x_sort(i+1))
countx=countx+1;
end
if (y_sort(i)~=y_sort(i+1))
county= county+1;
end
if (z_sort(i)~=z_sort(i+1))
countz= countz+1;
end
end
countx=993
county=1894
countz=1900
cutlength=countx*county*countz
marglen=[countx county countz];
X_reshape = reshape(x(1:cutlength),marglen);
Y_reshape = reshape(y(1:cutlength),marglen);
Z_reshape = reshape(z(1:cutlength),marglen);
U_reshape = reshape(U(1:cutlength),marglen);
result = trapz( squeeze(X_reshape(:,1,1)), trapz(squeeze(Y_reshape(1,:,1)),trapz(squeeze(Z_reshape(1,1,:)), U_reshape, 3), 2), 1);
disp(result);
But this time it yields out of memory error when i went for reshaping.
Can you help?
Walter Roberson

Walter Roberson (view profile)

on 19 May 2019
Have you considered
countx = length(unique(x));
county = length(unique(y));
countz = length(unique(z));
You overwrite your countx, county, countz. The values you assign multiply to 3573409800 . For the vector 1 to that, allow 8 bytes per entry, giving 28587278400 bytes, which is about 27 gigabytes. Writing the code that way also implies that x and y and z and U are also all at least that large, getting you up to over 108 gigabytes.