Extract data from countour and add two (or more) contours together.

My problem is thus:

I have a mesh onto which I 'plot' a function.

That is, for every point x,y, I have a value Q

What I want to do is extract the value of Q for every x,y and ADD (superimpose) it to the value of Q for every x,y .

My code plots a stress field centered at at point (0,0) then it is moved a length dx and dy and an appropratiate step length. Then this is plotted.

If I have more than one point, I wish to superimpose the data.

i.e. I wish to set up a loop where the original position is kept, then another stress field is generated at different dx, dy, and then the two are added.

regards William

u = 80000; % Shear modulus, measure in Mpa nu = 0.3; % Poissions ratio of material b = .00025; % Burgers Vector measured in microns

[x,y]=meshgrid(-100:1:100,-100:1:100); % area of grid to which stress field is mapped
angle = 54.7; % angle of slip plane in degrees
theta = angle*pi/180; % angle converted to radians
steplength = 400;
interplanarspacing = 100*b; % spacing between slip planes is about 100 b

x2 = (x.*cos(theta) + y.*sin(theta));% rotates and displaces dislocation
y2 = (y.*cos(theta) - x.*sin(theta));%; rotates and displaces dislocation

QxxPRIME = -(u .* b / (2.*pi .* (1-nu))) .* ((y2 .* (3.*x2.^2 +y2.^2)) ./ ((x2.^2 + y2.^2).^2)); % Qyy stress PRIME
QyyPRIME = (u .* b / (2.*pi .* (1-nu))) .* ((y2 .* (x2.^2 -y2.^2)) ./ ((x2.^2 + y2.^2).^2)); % Qyy stress PRIME
QxyPRIME = (u .* b / (2.*pi .* (1-nu))) .* ((x2 .* (x2.^2 -y2.^2)) ./ ((x2.^2 + y2.^2).^2)); % Qxy stress PRIME

QxxPRIME(QxxPRIME<-.5) = -.5; % erase low value
QxxPRIME(QxxPRIME>.5) = .5; % erase high value
dx = steplength.*cos(theta);
dy = steplength.*sin(theta);
y = y+dy;
x = x+dx;

caxis([-0.5,0.5]); % set the legend limit
contourf(x(1,:),y(:,1)',QxxPRIME)
axis equal

Answer 1

[Thorsten]

The idea is to compute a new QxxPRIME for the different offsets and add them to QxxPRIME0 for zero offset:

 u = 80000; % Shear modulus, measure in Mpa 
 nu = 0.3; % Poissions ratio of material 
 b = .00025; % Burgers Vector measured in microns

 N = 500; % grid size 

 [x,y]=meshgrid(-N:1:N,-N:1:N); % area of grid to which stress field is mapped
angle = 54.7; % angle of slip plane in degrees
theta = angle*pi/180; % angle converted to radians
steplength = 400;
interplanarspacing = 100*b; % spacing between slip planes is about 100 b
Nsteps = 10;

c = -(u*b/(2*pi*(1 - nu))); % constants independent of loop variable sl

for sl = linspace(0, steplength, Nsteps)

     dx = sl.*cos(theta);
     dy = sl.*sin(theta);
     x2 = x.*cos(theta) + y.*sin(theta) + dx;% rotates and displaces dislocation
     y2 = y.*cos(theta) - x.*sin(theta) + dy;%; rotates and displaces dislocation

     QxxPRIME = c*(y2.*(3*x2.^2 + y2.^2))./((x2.^2 + y2.^2).^2); % Qyy stress PRIME
     QxxPRIME(QxxPRIME<-.5) = -.5; % erase low value
     QxxPRIME(QxxPRIME>.5) = .5; % erase high value

     subplot(1,2,1)
     contourf(x,y,QxxPRIME)
     caxis([-0.5,0.5])
     axis square
     title(['Contour at Steplength ' num2str(sl)])

     if sl == 0
       QxxPRIME0 = QxxPRIME;
     else
       QxxPRIME = QxxPRIME0 + QxxPRIME;
     end

      subplot(1,2,2)
     contourf(x,y,QxxPRIME)
     caxis([-0.5,0.5])
     axis square
     title(['Steplength 0 + ' num2str(sl)])
     pause
  end

[William White]

Hi Thanks for the response, but this is not quite what I was after.

This has plotted two separate fields on the same graph; i.e. the values of Q at a first set of (x,y) are not added to the second values of Q and the second set of (x,y).

For example, Say that the value of Q1 at x=100, y=100 is 100 I want this value saved so that if The value of Q2 at x=100, y=100 is -30 Then this value is saved and Qtotal is Q1+Q2 so that the value is 70

And this is the case for all values of Q at all x,y.

hope this makes sense

[William White]

Hi Again I really do appreciate this, I am not sure this will work (?)

The formula for Q only works when the origin x,y = 0. i.e. the value of Q is infinite at zero

If the origin of Q is non-zero, then the value of Q at the orgin is non-zero and all the other values of Q at all x,y are 'out' by the same margin.

For example Imagine two fields generated by the formula for Q

One is centered at 10,10 The other is centered at 20,20

(both origins taken at random)

The stress fields defined by Q for both will be exactly the same, however, the interest is when you add them.

For example the values for Qtotal at 15,15 (for arguments sake) will be the values of each Q at that point, which if the function was symmetrical would be like the value for Q1 at 5,5 and the values of Q2 at -5,-5 added together.

Effectively, I am trying to superimpose the output of lots of Qs all at different locations. However, the field generated by Q1 at any x,y, has to be exactly the same as the field generated by Q2 at any other x,y

I hope this makes sense, I really do appreciate your input.

I did try to run the code, and got an error

Undefined function or variable 'QxxPRIME0'.

Error in test (line 28) QxxPRIME = QxxPRIME0 + QxxPRIME;

Sorry, I am very new to matlab and having to learn it quickly. Things that seem to be obvious to work out on paper are fiendish in code,

regards William

[Thorsten]

Hi William, I have updated the code. There is now a plot of the new displaced QxxPRIME that looks exactly as the QxxPrime at (0,0), and a plot of the addition. The plots are shown on a larger grid of size -500:500 such that the displaced plot is visible.

I ran the posted code without any error. The variable QxxPRIME0 is set at the first iteration of the loop, when sl takes the value 0. I do not understand how you can get the "Undefined function or variable 'QxxPRIME'." error.

[William White]

That looks very promising. Do you mind explaining through what is happening - particularly, how you got each stress field to be accurate at a non-zero origin?

Ultimately, my aim to to have lots (potentially thousands) of dislocations (the stress field is generated by a dislocation) in the plane, and to calculate the total stress field as the dislocations move.

For example, in the first instance, is it possible to set up a stress field centred at, say 10,10 and another at 15,15 and plot them on the same contour map? The positions of the dislocations is utterly important. One can lie on a line at 0°, or 54.7° or - 54.7 degrees but it can ONLY move along that lines. So if there are two on parallel lines at 54.7°, they can move towards each other, then apart, but never touch. They can only touch if they are on the same line.

regards William

[Thorsten]

The field is the same because the arguments x and y are just rotated and displaced, i.e., moved in the plane. If a function(x,y) is zero for (0,0), then the function f(x+dx, y+dy) is zero at x = -dx, y = -dy.

[William White]

so would it be possible to just plot two fucntions; one centred at 10,10 and another at say, 12,10? So that indvidually, they look exactly the same, but superimposed, give a new field?

Also, I think that by adding dx to the rotation, it has broken the correct working of the code?

For example, if my steplength is, say, 10; then for the stress field to move correctly, it moves 10 units along a 54.7 incline. So dx is 5.7 and dy is 8.2; that is, if the steplength is zero, then the origin of the stress field is at 0,0; if the steplength is 10, the new orign is at 5.7, 8.2.

I took the liberty of removing the loop from your code to see what was going on with the movement.

kind regards

u = 80000; % Shear modulus, measure in Mpa 
 nu = 0.3; % Poissions ratio of material 
 b = .00025; % Burgers Vector measured in microns

   N = 50; % grid size 

[x,y]=meshgrid(-N:.1:N,-N:.1:N); % area of grid to which stress field is mapped
angle = 54.7; % angle of slip plane in degrees
theta = angle*pi/180; % angle converted to radians
steplength = -10;
interplanarspacing = 100*b; % spacing between slip planes is about 100 b
Nsteps = 30;

c = -(u*b/(2*pi*(1 - nu))); % constants independent of loop variable sl

       dx = steplength*cos(theta);
       dy = steplength*sin(theta);
       x2 = (x.*cos(theta) + y.*sin(theta)) +dx;% rotates and displaces dislocation
       y2 = (y.*cos(theta) - x.*sin(theta)) +dy;%; rotates and displaces dislocation

       QxxPRIME = c*(y2.*(3*x2.^2 + y2.^2))./((x2.^2 + y2.^2).^2); % Qyy stress PRIME
       QxxPRIME(QxxPRIME<-.5) = -.5; % erase low value
       QxxPRIME(QxxPRIME>.5) = .5; % erase high value

       contourf(x,y,QxxPRIME)
       caxis([-0.5,0.5])
       axis equal
       title('Contour at Steplength ')

regards Willam