MATLAB Answers

Help indexing: How to index at the interfaces of cells and not the cells

1 view (last 30 days)
I am struggling to understand how to setup some indexing problem in MATLAB. If my entire domain is defined from 1 to N with a padding of two cells. So say my internal domain (domain without padded cells) from 2 to N-2, then if I want to calculate some values at the "interfaces" of the cells instead at the cells, how can I setup this on MATLAB? For example I have the following code for the internal domain that I would like to change my indicies to the interfaces?
My set up is like:
% calculate Flux at j+1/2
%
% j+1/2 Cell's grid:
% | wL| |
% | /|wR | 1 2 3 4 N-2 N-1 N
% | / |\ | {x=0} |-o-|-o-|-o-|-o-| ... |-o-|-o-|-o-| {x=L}
% |/ | \ | 1 2 3 4 5 N-1 N N+1
% | | \|
% | | | NC: Here cells 1 and N are ghost cells%
% j j+1
I keep getting the error:
Index in position 2 exceeds array bounds
The code for example:
function [ U_L, U_R] = MUSCL(U,N,NEQ )
size_U = size(U);
U_L = zeros(size_U);
U_R = zeros(size_U);
delta = 1e-6; % epsilon
for i= 1:NEQ
for j = 3:N-1
denom_L = (U(i,j) - U(i,j-1));
denom_R = (U(i,j+2) - U(i,j+1)); % error here since index exceeds the bounds of the array N (1001 here)
denom_L = sign(denom_L+1e-64)*max(delta, abs(denom_L));
denom_R = sign(denom_R+1e-64)*max(delta, abs(denom_R));
r_L = (U(i,j) - U(i,j-1))/denom_L;
r_R = (U(i,j+2) - U(i,j+1))/denom_R;
theta_L = Limiter(r_L);
theta_R = Limiter(r_R);
end
end
for j = 2:N-2
U_L(:,j) = U(:,j)+0.5*theta_R.*(U(:,j)-U(:, j-1));
U_R(:,j) = U(:,j+1)-0.5*theta_L.*(U(:,j+2)-U(:,j+1));
end
end
How can I find U_L, U_R, r_L, r_R at j+1/2 ??

Accepted Answer

per isakson
per isakson on 8 May 2021
Edited: per isakson on 8 May 2021
Index in position 2 exceeds array bounds
U(:,j+2) is most likely the problem.
I assume that the size of U is NEQxN. The maximum value of j is N-1 (for j = 3:N-1). Thus, the maximum value of j+2 is N+1, which exceeds the array bound, N
  2 Comments
per isakson
per isakson on 9 May 2021
Your code is about "flux", "cells" and "interfaces". The name, MUSCL, does that stand for "muscle"? However, I neither understand the physics nor the code. Futhermore, I cannot run the code, since I neither have the function, Limiter(), nor a realistic sample of U.
It looks strange to me that none of the left-hand-side variables in the double loop is indexed.
I've made a little demo in response to your answer to @DGM.
%% A sample matrix
M = 40 .* membrane( 1, 25 );
%%
to_the_interfaces = @(row) interp1( (1:numel(row)), row, (1.5:1:numel(row)-0.5) );
row = (1:10);
to_the_interfaces(row)
ans = 1×9
1.5000 2.5000 3.5000 4.5000 5.5000 6.5000 7.5000 8.5000 9.5000
%%
M_at_the_interfaces = nan(size(M)-[0,1]);
for jj = 1 : size(M,1)
M_at_the_interfaces(jj,:) = to_the_interfaces(M(jj,:));
end
imagesc( M_at_the_interfaces )

Sign in to comment.

More Answers (1)

DGM
DGM on 8 May 2021
Edited: DGM on 8 May 2021
Just interpolate.
x = 0:10;
y = x.^2
xhalf = 0.5:9.5;
yhalf = interp1(x,y,xhalf,'linear')
plot(x,y,'rd',xhalf,yhalf,'bo')
  1 Comment
Lujain Almarhabi
Lujain Almarhabi on 8 May 2021
I have no clue how to apply this in my example. Can you show something for my for loop?

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!

Translated by