I want write the following matrix as a sparse matrix. I want to run my code for M=N=K=100. However, my computer always ran out of memory the way I define the matrix.

1 view (last 30 days)
Kumu
Kumu on 11 Jul 2019
Edited: Hari Krishna Ravuri on 1 Aug 2019
M=8;
N=8;
K=8;
m=M-1;%number of interior grid points
n=N-1;%number of exterior grid points
k=K-1;%number of exterior grid points
%%%USER DEFINE CONSTANTS AND FUNCTIONS%%%%%
%define the rectangular domain
a=-1;
b=1;
c=-1;
d=1;
e=-1;
f=1;
%defining beta coefficient
betamin=1;% inside the interface
betaplus=1;%outside the interface
row=betamin/betaplus;
sigmaplus=0;
sigmamin=0;
utrue_m = @(x,y,z) exp(x.^2+y.^2+z.^2);
utrue_p = @(x,y,z) exp(x.^2+y.^2+z.^2)+1;%%%USER DEFINE CONSTANTS AND FUNCTIONS CLOSE%%%%%
phi =@(x,y,z) x*x+y*y+z*z-1/4; %interface is x*x+2*y*y+z*z-1/4;
h=(b-a)/M;
hs=h^2;
mnk=m*n*k;
mn=m*n;
x=a+(1:m)*h; % set mesh values in x-direction
y=c+(1:n)*h;
z=e+(1:k)*h;
n1 = 0;
index=zeros(m,n,k);% memory allcation
digits(32)
for i=2:m-1
for j=2:n-1
for l=2:k-1
index(i,j,l) = 0 ; % regular grid point
if(phi(x(i),y(j),z(l)) == 0)
index(i,j,l) = 1; % on the interface
end
if(phi(x(i),y(j),z(l)) > 0) %outside the interface
if(phi(x(i),y(j),z(l))*phi(x(i-1),y(j),z(l)) <= 0 || phi(x(i),y(j),z(l))*phi(x(i+1),y(j),z(l)) <= 0)%irregular point
index(i,j,l) = 1;
end
if(phi(x(i),y(j),z(l))*phi(x(i),y(j-1),z(l)) <= 0 || phi(x(i),y(j),z(l))*phi(x(i),y(j+1),z(l))<= 0) %irregular point
index(i,j,l) = 1;
end
if(phi(x(i),y(j),z(l))*phi(x(i),y(j),z(l-1)) <= 0 || phi(x(i),y(j),z(l))*phi(x(i),y(j),z(l+1))<= 0) %irregular point
index(i,j,l) = 1;
end
end
if(phi(x(i),y(j),z(l)) < 0) %inside the interface
if(phi(x(i),y(j),z(l))*phi(x(i-1),y(j),z(l)) < 0 || phi(x(i),y(j),z(l))*phi(x(i+1),y(j),z(l)) < 0)%irregular point
index(i,j,l) = 1;
end
if(phi(x(i),y(j),z(l))*phi(x(i),y(j-1),z(l)) < 0 || phi(x(i),y(j),z(l))*phi(x(i),y(j+1),z(l))< 0) %irregular point
index(i,j,l) = 1;
end
if(phi(x(i),y(j),z(l))*phi(x(i),y(j),z(l-1)) < 0 || phi(x(i),y(j),z(l))*phi(x(i),y(j),z(l+1))< 0) %irregular point
index(i,j,l) = 1;
end
end
if(index(i,j,l) == 1)
n1 = n1 + 1;
index(i,j,l) = n1;% indexing all irregular points
end
end
end
end
for i = 1:m
for j = 1:n
for l = 1:k
if index(i,j,l)<1
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+(l-1)*mn) = -6 ;
if i~=m
A(i+(j-1)*m+(l-1)*mn,(i+1)+(j-1)*m+(l-1)*mn) = 1;
end
if i~=1
A(i+(j-1)*m+(l-1)*mn,(i-1)+(j-1)*m+(l-1)*mn) =1;
end
if j~=1
A(i+(j-1)*m+(l-1)*mn,i+((j-1)-1)*m+(l-1)*mn) = 1;
end
if j~=n
A(i+(j-1)*m+(l-1)*mn,i+((j+1)-1)*m+(l-1)*mn) = 1;
end
if l~=1
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+((l-1)-1)*mn) = 1;
end
if l~=k
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+((l+1)-1)*mn) = 1;
end
else
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+(l-1)*mn)= -6 +.14;
if i~=1 && j~=l && l~=l
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j-1)-1)*m+((l-1)-1)*mn) = 0.1;
end
if j~=1 && l~=1
A(i+(j-1)*m+(l-1)*mn,i+((j-1)-1)*m+((l-1)-1)*mn)= 0.2;
end
if i~=m && j~=1 && l~=1
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j-1)-1)*m+((l-1)-1)*mn)= 0.3;
end
if i~=1 && l~=1
A(i+(j-1)*m+(l-1)*mn,(i-1)+(j-1)*m+((l-1)-1)*mn)= 0.4;
end
if l~=1
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+((l-1)-1)*mn)= 1+5;
end
if i~=m && l~=1
A(i+(j-1)*m+(l-1)*mn,(i+1)+(j-1)*m+((l-1)-1)*mn)= 0.6;
end
if i~=1 && j~=n && l~=1
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j+1)-1)*m+((l-1)-1)*mn)= 0.7;
end
if j~=n && l~=1
A(i+(j-1)*m+(l-1)*mn,i+((j+1)-1)*m+((l-1)-1)*mn)= 0.8;
end
if i~=m && j~=n && l~=1
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j+1)-1)*m+((l-1)-1)*mn)= 0.9;
end
if i~=1 && j~=1
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j-1)-1)*m+(l-1)*mn)= 10;
end
if j~=1
A(i+(j-1)*m+(l-1)*mn,i+((j-1)-1)*m+(l-1)*mn)= 1+.11;
end
if i~=m && j~=1
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j-1)-1)*m+(l-1)*mn)= .12;
end
if i~=1
A(i+(j-1)*m+(l-1)*mn,(i-1)+(j-1)*m+(l-1)*mn)= 1+.13;
end
if i~=m
A(i+(j-1)*m+(l-1)*mn,(i+1)+(j-1)*m+(l-1)*mn)= 1+.15;
end
if i~=1 && j~=n
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j+1)-1)*m+(l-1)*mn)= .16;
end
if j~=n
A(i+(j-1)*m+(l-1)*mn,i+((j+1)-1)*m+(l-1)*mn)= 1+.17;
end
if i~=m && j~=n
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j+1)-1)*m+(l-1)*mn)= .18;
end
if i~=1 && j~=1 && l~=k
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j-1)-1)*m+((l+1)-1)*mn)= .19;
end
if j~=1 && l~=k
A(i+(j-1)*m+(l-1)*mn,i+((j-1)-1)*m+((l+1)-1)*mn)= .20;
end
if i~=m && j~=1 && l~=k
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j-1)-1)*m+((l+1)-1)*mn)= .21;
end
if i~=1 && l~=k
A(i+(j-1)*m+(l-1)*mn,(i-1)+(j-1)*m+((l+1)-1)*mn)= .22;
end
if l~=k
A(i+(j-1)*m+(l-1)*mn,i+(j-1)*m+((l+1)-1)*mn)=1+.23;
end
if i~=m && l~=k
A(i+(j-1)*m+(l-1)*mn,(i+1)+(j-1)*m+((l+1)-1)*mn)= .24;
end
if i~=1 && j~=n && l~=k
A(i+(j-1)*m+(l-1)*mn,(i-1)+((j+1)-1)*m+((l+1)-1)*mn)= .25;
end
if j~=n && l~=k
A(i+(j-1)*m+(l-1)*mn,i+((j+1)-1)*m+((l+1)-1)*mn)= .26;
end
if i~=m && j~=n && l~=k
A(i+(j-1)*m+(l-1)*mn,(i+1)+((j+1)-1)*m+((l+1)-1)*mn)=.27;
end
end
end
end
end

Sign in to answer this question.

Answers (1)

Hari Krishna Ravuri
Hari Krishna Ravuri on 24 Jul 2019
Edited: Hari Krishna Ravuri on 1 Aug 2019
I see that you are modifying the array A nested in three for loops and you are planning to run the entire script with M=100 , N=100 , K=100.With these values, the instructions in the three nested for loops will run for exactly 100x100x100 times which takes considerable amount of CPU time. For every iteration, the array A is being modified, as the Array was not initialized before, this will involve re-allocating the memory for the array and copying the data from the previously allocated array. Array of that size requires memory of approximately 7014GB which exceeds maximum array size preference. Creation of arrays greater than this limit may take a long time and cause MATLAB to become unresponsive. Please refer https://www.mathworks.com/help/matlab/matlab_env/set-workspace-and-variable-preferences.html for more information. It's evident that your algorithm has space complexity of O(n^3) which is very costly. You can optimize your algorithm to reduce the space complexity.

Categories

Find more on Creating and Concatenating Matrices in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!