Effacer les filtres
Effacer les filtres

for loop and function handle

11 vues (au cours des 30 derniers jours)
Jorge Garcia Garcia
Jorge Garcia Garcia le 13 Juil 2023
Hello, I would really appreciate an answer as this is driving me crazy.
I have a system of PDE, and I have a plate geometry. I have stored the values of the edges in the perimeter in one vector and the values of the inner edges in another vector.
When I apply boundary conditions to the perimeter:
in=0;
for k = 1:numel(perimeter)
edge = perimeter(k);
in=in+1;
bcfunc= @(location,state) mybc(location,state,dis,ts,in);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",[edge,edge],"u",bcfunc,"EquationIndex",1,"q",[0 0],"g",0);
end
in=0;
where
function bcMatrix = mybc(location,state,acs,ti,i)
T=state.time;
% Check if T is NaN and assign the previous value if true
if isnan(T)
%T = previousTime;
vq = NaN(size(location.x));
else
%previousTime = T; % Update previousTime with the current value of T
ac = acs(i,:);
vq = interp1(ti, ac, T);
end
bcMatrix = vq;
end
it works. I pass dis, a matrix containing rows that represent the displacement of every edge in the perimeter, ts the time measurement and in and index to move along the rows depending on the edge.
Now instead of the perimeter I want to do the same in the inner edges, the values stored in inedges
in=0;
for k = 1:numel(inedges)
edge = inedges(k);
in=in+1;
bcfunc= @(location,state) mybc(location,state,dis,ts,in);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",[edge,edge],"u",bcfunc,"EquationIndex",1,"q",[0 0],"g",0);
end
in=0;
and in this ocasion, the corresponding edges do not enter into mybc, so it does not apply the boundary conditions. Am I missing something here?This is driving me crazy as it does not make sense

Réponse acceptée

ProblemSolver
ProblemSolver le 13 Juil 2023
To my understanding of the code snippet that you provided it can be couple of things. To begin with you can try these things:
  • Verify if the 'inedges' vector is correct. You can use this before each loop.
disp(inedges); % Print the contents of the inedges vector
  • Print a statement of just show if the loop is correct or not something like this. In this way you which edge is being currntly processed and if the entry is correct or not.
for k = 1:numel(inedges)
edge = inedges(k);
disp(edge); % Print the current edge value
in = in + 1;
bcfunc = @(location,state) mybc(location,state,dis,ts,in);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",[edge,edge],"u",bcfunc,"EquationIndex",1,"q",[0 0],"g",0);
end
  • You can either use MATLAB debugger or this snipped provided below to check if your function is providing the desired answers or not. If not, then you know where you are making the mistake.
function bcMatrix = mybc(location,state,acs,ti,i)
T = state.time;
% Print the inputs for debugging
disp(location);
disp(state);
disp(acs);
disp(ti);
disp(i);
% Rest of your code...
bcMatrix = vq;
end
Based on the details you provided, I could only help you so much. I hope these suggestions help you to figure out where you are going wrong. You should also read this debugging documentation, so that you can apply them as required:
  4 commentaires
ProblemSolver
ProblemSolver le 14 Juil 2023
In you updated code:
You are trying to pass the variables dis1 and ts to the function mybc using the anonymous function. Please read this link to understand the definition of function in MATLAB (https://www.mathworks.com/help/matlab/ref/function.html). You have send three input points 'location', 'state', and 'dis1'. Therefore modify the function accordingly.
Secondly, defining the anonymous function (https://www.mathworks.com/help/matlab/matlab_prog/anonymous-functions.html). PLease read and understand it. You have to define the inuts arguments correctly. The 'mybc' function expects four input arguments, therefore include them as well.
Please post your code properly and questions and explanation on what you have done. Read this: https://www.mathworks.com/matlabcentral/answers/6200-tutorial-how-to-ask-a-question-on-answers-and-get-a-fast-answer
Jorge Garcia Garcia
Jorge Garcia Garcia le 17 Juil 2023
sorry, it is that so many trials to see what may be wrong, made me copy the wrong code. here is is the whole code to see whether it helps to understand:
%-----------Define dimensions of the plate-------------------------------
len = 3; % length in x direction
width = 2; % Width in y direction
nodes = 1.22; % Distance between the different horizontal and vertical battens
edges=[];
% Create the geometry description matrix-----------------------------------
region_boundaries_y = [0];
dist_y = 0;
interval = nodes;
while dist_y < width
dist_y = dist_y + interval;
if dist_y > width
dist_y = dist_y - (dist_y - width);
elseif abs(dist_y - width) < 1e-10
dist_y = width;
end
region_boundaries_y = [region_boundaries_y dist_y];
end
% Create divisions in X (horizontal divisions) ----------------------------
region_boundaries_x = [0];
dist_x = 0;
interval = nodes;
while dist_x < len
dist_x = dist_x + interval;
if dist_x > len
dist_x = dist_x - (dist_x - len);
elseif abs(dist_x - len) < 1e-10
dist_x = len;
end
region_boundaries_x = [region_boundaries_x dist_x];
end
num_regions_x = length(region_boundaries_x) - 1;
num_regions_y = length(region_boundaries_y) - 1;
%---------Create the PDE Model with a single domain----------------------
model=createpde;
%-CREATE THE GEOMETRY MATRIX---------------------------------------------
% Coordinates of the boundaries in x and y--------------------------------
coords_x = region_boundaries_x';
coords_y = region_boundaries_y';
sq=0.05;
% Create the geometry description matrix
gdm = [];
for k = 1:numel(coords_x)-1
for j = 1:numel(coords_y)-1
% Check if the coordinates are completely inside the external plate
if coords_x(k) >= 0 && coords_x(k) <= len && coords_y(j) >= 0 && coords_y(j) <= width
% Add the vertical dimension
gdm = [gdm; 2, 5, coords_x(k), coords_x(k+1), coords_x(k+1), coords_x(k+1), coords_x(k), coords_y(j), coords_y(j), coords_y(j)+sq, coords_y(j+1), coords_y(j+1)];
end
end
end
%--------------------------------------------------------------------------
gdm=gdm';
%---------------Create the geometry--------------------------------------
[dl,bt] = decsg(gdm);
for j = 1:numel(region_boundaries_x)
for i=1:numel(region_boundaries_y)
indices = find(dl(2,:) == region_boundaries_x(j) & dl(3,:) == region_boundaries_x(j)&...
dl(4,:) == region_boundaries_y(i) & dl(5,:) == region_boundaries_y(i)+sq);
edges = [edges, indices];
end
end
perimeter = find((dl(2,:)== 0) & (dl(3,:) == 0) |...
(dl(4,:) == 0) & (dl(5,:) == 0)|...;
(dl(4,:) == width) & (dl(5,:) == width)|...
(dl(2,:)== len) & (dl(3,:) == len));
%Allows to select just the left edges of the inner squares
edges = setdiff(edges, perimeter);
geometryFromEdges(model,dl);
% Generate the mesh and plot the geometry---------------------------------
msh = generateMesh(model);
%---------Create the different regions of the model-----------------------
% Obtain nodes and elements of the mesh
[p,e,t] = meshToPet(model.Mesh);
modelTwoDomain = createpde(2);
geometryFromEdges(modelTwoDomain,dl);
%
% Now we create the PDE systems as symbolic equations---------------------
syms pres
syms u1(x,y,t) u2(x,y,t)
pdeeq = [-laplacian(u1,[x y])+u2; D*laplacian(u2,[x y])+ mass*diff(u1,t,t)-pres];
symcoeffs = pdeCoefficients(pdeeq,[u1,u2],'Symbolic',true);
c2=symcoeffs.c;
m2=symcoeffs.m;
symcoeffs=subs(symcoeffs,pres,1);
coeffs=pdeCoefficientsToDouble(symcoeffs);
specifyCoefficients(modelTwoDomain,'m',coeffs.m,'d',coeffs.d,'c',...
coeffs.c,'a',coeffs.a,'f',coeffs.f);
% INITIAL CONDITIONS ---------------------------------------------------
setInitialConditions(modelTwoDomain,[0;0],[0;0]);
dis=12x13000 it is this is the dimension of the matrix containing displacements. Basically I want that the first edge in the loop gets the first row of the matrix, and so on.
tim=[ts(20) ts(30) ts(35)];
for k = 1:numel(edges)
edge = edges(k);
dis1=dis(k,:);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",edge,"u",@(location,state)mybc1(location,state,dis1,ts),"EquationIndex",1);
end
res=solvepde(modelTwoDomain,tim);
% Access the solution at the nodal locations
sol=res.NodalSolution;
In bold is the part of the code that gives me problems. the vector perimeter contains the values: 1 2 7 8 9 13 14 15 16 17 18 19
if I use
for k = 1:numel(perimeter)
edge = perimeter(k);
dis1=dis(k,:);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",edge,"u",@(location,state)mybc1(location,state,dis1,ts),"EquationIndex",1);
end
where
function bcMatrix = mybc1(location,state,acs,ti)
T=state.time;
if isnan(T)
%T = previousTime;
vq = NaN(size(location.x));
else
ac = acs;
vq = interp1(ti, ac, T);
end
bcMatrix = vq;
end
The code runs without not problem. U stablish a breaking point in the function in space.time, and the code runs into it.
On the other hand if I use what I need:
for k = 1:numel(edges)
edge = edges(k);
dis1=dis(k,:);
applyBoundaryCondition(modelTwoDomain,"mixed","Edge",edge,"u",@(location,state)mybc1(location,state,dis1,ts),"EquationIndex",1);
end
then the code does not enter the function. Vector edges contains the edges I need edges = 3 5 20 22
I really cannot see why in one case it is entering the function mybc1 and in the other case of edges, it is not.

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Data Distribution Plots dans Help Center et File Exchange

Produits


Version

R2023a

Community Treasure Hunt

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

Start Hunting!

Translated by