I am working with a PDE that has one independent variable and varies only in space(m = 0, d = 0, only one equation). I have a 2D domain with multiple subdomains. Every subdomain has its own c, a, and f coefficients. At first, I used the specifyCoefficients() function in a loop over ever face/subdomain, assigning to each one its own constant coefficients. Now I am attempting to use a function to define those same coefficients, but the number of regions reported by the coefficients (using disp(length(region.x)) inside each function) is not consistent between c, f, and a, and is not consistent with the number suggested by the MathWorks website (see below). The contour plots of the two different solutions are attached (the only difference between them is the method I used to assign "a"), and below I will list the specifics of my pde object, as well as relevant code snippets.
Number of faces = 200
Number of mesh nodes = 4890
Number of elements = 9532 ('linear' mesh)
Code used when assigning coefficients:
[geom_setup,geom_bt] = decsg(gdm,g_formula,g_name);
geom_return = geometryFromEdges(Model,geom_setup);
[point_info,edge_info,tri_info] = meshToPet(Model.Mesh);
f_sub_num = tri_info(end,:)';
f_elements = length(f_sub_num);
for i = 1:f_elements
face_num = f_sub_num(i);
voxel_num = find(geom_bt(face_num,:),1);
ii = g_pts_y - (mod(voxel_num,g_pts_y) +...
g_pts_y*(mod(voxel_num,g_pts_y) == 0)) + 1;
jj = floor((voxel_num-1)/g_pts_y) + 1;
CM_col = CM(:,data(ii,jj)==CM(1,:));
f_sub_coef(i) = CM_col(4);
for i = 1:g_pts_y
for j = 1:g_pts_x
voxel_num = (g_pts_y-(i-1)) + g_pts_y*(j-1);
CM_col = CM(:,data(i,j)==CM(1,:));
c = CM_col(2);
a = CM_col(3);
f = @F_coeffunction;
face_num = find(geom_bt(:,voxel_num),1);
function [f] = F_coeffunction(region,state)
disp(['F length ' num2str(length(region.subdomain))]);
f = f_sub_coef(1:length(region.subdomain));
Note that length(region.x) = 9532 (the number of elements). However, if I attempt to use the same code seen in F_coefficient, for the "a" coefficient (the only difference is the value read from the CM_col plus an extra global variable so that "f" and "a" do not overwrite each other's values), I run into the following problem. When length(region.x) is called inside of a, it returns 3x9532, which corresponds to numel(Model.Mesh.Elements), not length(Model.Mesh.Elements(1,:))! The website's instructions say "Nr is equal to the length of the region.x or any other region field" for ALL coefficients a, c, and f.
Also, while running solvepde(Model), I noticed that length(region.x) changes from 1 to 9532. I thought that region.x was strictly related to the number of elements in the mesh, which caused problems when I tried to assign coefficients to each triangle as shown above. This led to the patching the code using f = "f_sub_coef(1:length(region.subdomain));".
QUESTION 1: Why does region.x change size for the f coefficient during the execution of solvepde(Model), and given that variation, how do I assign f coefficients for each triangle?
QUESTION 2: Why does region.x change size for "a" and "f"?
Secondly, my goal is to assign different a, c, and f values for each subdomain. It is easy to go back and forth between mesh triangles and face numbers using the outputs of meshToPet(). However, I have tried creating a coefficient vector for a using the code,
a_sub_coef = zeros(3,f_elements);
a_sub_coef = reshape(a_sub_coef,1,numel(a_sub_coef));
This way, there is a 1-1 correspondence between the columns of a_sub_coef and the triangles in the Mesh. I then assume that, after reusing reshape(), there is a 1-1 correspondence between the elements of a_sub_coef and region.x. However, when executing this code I obtain results that are different from the first approach.
QUESTION 3: How should I write @A_coefficient(region,state) so that I perform an assignment in the exact same fashion as looping over the faces and using specifyCoefficients(...)?
Recall, each face (therefore, the mesh triangles corresponding to those faces) should each have their own constant values of a, c, and f.
Thank you for your help with this very involved question!!