Store tensor history into multidimensional array

7 views (last 30 days)
Dear all,
I read the documentatuin on Multidimensional Arrays and I haven't been able to find a solution to my problem.
I have a tensor history like this:
s11 = [xx1, xx2, xx3];
s22 = [yy1, yy2, yy3];
s33 = [zz1, zz2, zz3];
s12 = [xy1, xy2, xy3];
s13 = [xz1, xz2, xz3];
s23 = [yz1, yz2, yz3];
I want to store these tensor component histories in a single multidimensional array TENSOR, such that each slice of the array is a complete stress tensor at a certain point in time, e.g.;
>>TENSOR(:, :, 1.0)
xx1 xy1 xz1
xy1 yy1 yz1
xz1 yz1 zz1
>>TENSOR(:, :, 2.0)
xx2 xy2 xz2
xy2 yy2 yz2
xz2 yz2 zz2
(and so on)
I know how to create a multidimensional array in individual slices. For example, I could simply define the first slice then use
TENSOR(:, :, n) = [s11(n), s12(n), s13(n); s12(n), s22(n), s23(n); s13(n), s23(n), s33(n)];
To add the next slide. However, I would like to define the array in a single assignment, since the whole point of the exercise is to avoid inefficient loops.
I would greatly appreciate any advice!
Best regards,
  1 Comment
Louis Vallance
Louis Vallance on 4 Dec 2016
Andrei's solution below works, although I had to make a couple of tweaks to get the behaviour just right so that it works in general for any 3x3xL multidimensional array:
Construct a 3x3xL multidimensional array of the
stress tensor history, where L is the history
% Direct stress components
normals = [Sxx; Syy; Szz]';
% Shear stress components
shears = [Txy; Txz; Tyz]';
% Assign normal stress component histories to diagonal terms
diagonals = repmat(eye([3.0, 3.0]), 1.0, 1.0, L);
diagonals(diagonals > 0.0) = normals';
% Assign shear stress component histories to non-diagonal terms
nonDiagonals = repmat(tril(ones(3.0), -1.0), 1.0, 1.0, L);
nonDiagonals(nonDiagonals > 0.0) = shears';
% Combine diagonal and non-diagonal terms
multidimensionalStressTensor = diagonals + nonDiagonals +
permute(nonDiagonals, [2.0, 1.0, 3.0]);
Get the Eigenvalues for each frame of the tensor data. This code
was written by Bruno Luong and can be found on the exchange (File
ID: #27680)
eigenvalues = eig3(multidimensionalStressTensor);
The principal stresses are the maximum, median and
minimum values of the Eigenvector at each tensor
s1 = max(eigenvalues);
s2 = median(eigenvalues);
s3 = min(eigenvalues);
Initial tests indicate that this method is approximately 4 times faster than a FOR loop with 2D matrices.

Sign in to comment.

Accepted Answer

Andrei Bobrov
Andrei Bobrov on 4 Dec 2016
sxx = reshape(1:9,3,[]); % sxx = [s11;s22;s33];
sxyz = reshape(11:19,3,[]); % sxyz = [s12;s13;s23]
s = size(sxx);
D = repmat(eye(s),1,1,s(2));
D(D>0) = sxx;
ND = repmat(tril(ones(s),-1),1,1,s(2));
ND(ND>0) = sxyz;
tenzor1 = D + ND + permute(ND,[2,1,3]);
  1 Comment
Louis Vallance
Louis Vallance on 4 Dec 2016
Thanks Andrei. This works as long as I make the following change to your code:
D(D>0) = sxx';
ND(ND>0) = sxyz';
Otherwise you get tensor histories being assigned to incorrect tensor locations. With the transpose modification it achieves exactly what I wanted. Thanks!

Sign in to comment.

More Answers (1)

Daniel kiracofe
Daniel kiracofe on 4 Dec 2016
As far as I am aware, it is not possible to do what you want in matlab. There is no such thing as a 3-dimensional array literal, which is what you want. It is possible in some programming languages, such as C++, but not matlab. The only way I can think of do it without loops is to create a 1D array and then reshape it. e.g:
X = [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27];
y = reshape(X, [3 3 3])
which would give you this output:
y =
ans(:,:,1) =
1 4 7
2 5 8
3 6 9
ans(:,:,2) =
10 13 16
11 14 17
12 15 18
ans(:,:,3) =
19 22 25
20 23 26
21 24 27
Depending on what you are trying to do, you may find this to be useful. I've not used it, and the documentation on their website is sparse, but it is a toolbox for tensor manipulation in matlab which may helpful
  1 Comment
Louis Vallance
Louis Vallance on 4 Dec 2016
I Daniel,
I think it's possible with the solution from Andrei. I'll check out the TensorToolbox though, because it sounds like something which would be very useful to me.
Thanks for the help!

Sign in to comment.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!