How to correctly implement a formula that describes the deformation of a plate

Hello,

I tried to implement the following formula into Matlab:

First I tried to do it with symbolic variables like this (obviously i put values on the parameters and declared the symbolic variables n, N, m, M, x, y and w):

w = ((16.*q)./(K.*pi.^6)).*symsum(subs(symsum(subs((1./(m.*n.*(((m.^2)./(a.^2)) + ((n.^2)./(b.^2))))).*sin(x.*sym(3.14159)./a).*sin(y.*n.*sym(3.14159)./b), n, 2.*N-1), N, 1, 10), m, 2.*M-1), M, 1, 10);

but when I tried to display w, it seemed that Matlab didn't calculate it all the way.

So i tried it the numerical way (again the parameters were declared before):

k = 41;
x1 = linspace(0, a, k);
x2 = linspace(0, b, k);
[X, Y] = meshgrid(x, y);

w = zeros(k,k);
for m=1:2:21
    for n=1:2:21
        w = w + 1./(m.*n.*(((m.^2)./(a.^2))+((n.^2)./(b.^2)))).*sin(X.*m.*pi./a).*sin(Y.*n.*pi./b) ;
    end
end
w = ((16*q)/(K*pi^6))*w;

I get results that way but they vary when I change the value for k. But k should only help with the meshgrid and the size of the matrix of w. Why does it change the values for w at the same place (for example the center)? What have I done wrong?

[Simson Hutagalung]

How to generate deformation colormap 2D in matlab with data from excel?

Answer 1

[Walter Roberson]

I get results that way but they vary when I change the value for k. But k should only help with the meshgrid and the size of the matrix of w. Why does it change the values for w at the same place (for example the center)

x1 = linspace(0, a, k);
x2 = linspace(0, b, k);
[X, Y] = meshgrid(x, y);

You do not define x or y for us, but you also do not use x1 or x2 later in your code, so we assume that x is one of x1 or x2 and that y is the other one.

sin(X.*m.*pi./a)

X is 0 to a, so X/a is a fraction 0 to 1. You multiply that by π so by itself that is defining a half trig cycle. When you then multiply by m, you are effectively creating a sine wave with frequency m/2

Now, suppose your linspace(0,a,k) can include "exactly" X = a/2 and then with the division by a that gives and m is an odd integer. That is (m-1)/2 full cycles and a left over half cycle. sin() of a half cycle is, to within round-off error, +/- 1. So if the linspace() can include "exactly" half way through a then the term will give a min or max

When does linspace(0,a,k) include a/2 ? Answer: when k is odd. The formula is (0:k-1)/(k-1) * (a-0) + 0 and if k is even then k-1 is even and (k-1)/2 is an integer, and divide it by (k-1) you get 1/2 to within round-off, so you would generate a/2 to within round-off error.

What happens if k is even? Then you do not generate the mid-point. For example k = 4, (0:4-1)/(4-1) -> (0:3)/3 -> 0/3, 1/3, 2/3, 3/3 and none of those multiplied by a can be a/2 . But if you do not generate the mid-point then the sin(X.*m.*pi./a) is not going to be taking sin() of an odd multiple of and you are not going to get that clean +/- 1.

So... changing k is expected to change the results.

[NL]

Thank you, I understand it better now. I'm sorry for the x1 and x2 (I used other names for the variables before) but your assumption is correct.

[Simson Hutagalung]

How to generate deformation colormap 2D in matlab with data from excel?

[Walter Roberson]

@Simson Hutagalung

That is not related to the current Question. You should ask a new Question, and give more detail. Ask

Answer 2

[David Hill]

a=10;b=15;q=5;k=41;%do not know values
[m,n]=meshgrid(1:2:1000);
w=@(x,y)16*q/k/(pi^6)*sum(1./(m.*n.*(m.^2/a^2 + n.^2/b^2).^2).*sin(pi*m*x/a).*sin(pi*n*y/b),'all');
[X, Y] = meshgrid(linspace(0, a, k), linspace(0, b, k));
W=zeros(k);
for ii=1:numel(X)
  W(ii)=w(X(ii),Y(ii));
end