2

I have 40 structures in my Workspace. I Need to write a script to calculate the directional derivatives of all the elements. Here is the code :

[dx,dy] = gradient(Structure_element_1.value);  
dxlb = min(min(dx));  
dxub = max(max(dx));  
dylb = min(min(dy));  
dyub = max(max(dy));  

[ddx,ddy] = gradient(gradient(Structure_element_1.value));
ddxlb = min(min(ddx));  
ddxub = max(max(ddx));
ddylb = min(min(ddy));
ddyub = max(max(ddy));  

This is the code for one element. I Need to find out the same for all the 40 elements and then use it later. Can anyone help with this.

Aakash Sehgal
  • 171
  • 1
  • 2
  • 12
  • 2
    Do the structures all contain the same fields? If you have `Structure_element_1, Structure_element_2, ...` then it would be much better to have `Structure_element(1), Structure_element(2), ...` which allows you to more easily operate over the entire range of your data. The former essentially requires the use of `eval`, which is not good practice for many reasons, chief among them that it cannot be optimized by MATLAB's compiler. – sco1 Feb 11 '16 at 15:15
  • You probably want to read up on `for` loops combined with `eval` (see here: http://www.mathworks.com/matlabcentral/newsreader/view_thread/275986). Also, do `min(dx(:))` instead of `min(min(dx))`. – Mad Physicist Feb 11 '16 at 15:16
  • Also, you can probably reuse the result from the first call to `gradient` to pass into the second call instead of doing `gradient(gradient(...));`. – Mad Physicist Feb 11 '16 at 15:18
  • @MadPhysicist the only reading needed on `eval` is [this answer](http://stackoverflow.com/questions/32467029/how-to-put-these-images-together/32467170#32467170) which explains rather adequately why you should never use it if possible. – Adriaan Feb 11 '16 at 15:24
  • @Adriaan. Your answer is definitely true in general. If the OP is doing a throwaway script, they may disregard your advice though. Using `eval` almost always means that you are doing something wrong, whether it's in MATLAB or Python or any equivalent. – Mad Physicist Feb 11 '16 at 16:21
  • 3
    @MadPhysicist the problem with encouraging people to use `eval` in "throwaway scripts" is that they will get used to naming variables like this and will keep on using `eval`, even in non-throwaway scripts. Besides, either using `eval` which basically breaks MATLAB, or copying the same lines loads and loads of times. Using storage containers like cell arrays or proper structures is easier, less typing, less prone to errors and better efficiency wise. – Adriaan Feb 11 '16 at 16:27
  • Fair enough. The correct answer is "don't do it" in this case. Go back to @excaza 's comment and use that information instead. The stuff I said about `min`, `max` and `gradient` still applies though. – Mad Physicist Feb 11 '16 at 16:41
  • 1
    My advice is to go back to the script where you generated these 40 variables and put them into a sensible array so you don't have this problem. – beaker Feb 11 '16 at 17:53
  • @Adriaan thank you for your suggestions. – Aakash Sehgal Feb 12 '16 at 08:35
  • @excaza yes all the structures have the same field. Could you Elaborate a bit more on your solution. – Aakash Sehgal Feb 12 '16 at 08:36
  • @MadPhysicist Thank you for the advice. But i didnt really understand the use of 'eval' – Aakash Sehgal Feb 12 '16 at 08:36
  • @MadPhysicist when i try to reuse the value of the 1st Gradient to calculate the seond Gradient, the values are different in both the cases. '[dx,dy] = gradient(Structure_element_1.value); ' and then '[ddx,ddy] = grad([dx,dy]);' – Aakash Sehgal Feb 12 '16 at 08:55
  • @AakashSehgal. After looking at http://www.mathworks.com/help/matlab/ref/gradient.html, I am pretty sure that your way is actually wrong. `[p q] = gradient(gradient(blah))` will acutally set `p = d^2/dx^2 blah` and `q = d/dy(d/dx blah)`. I am pretty sure you want `q = d^2/dy^2 blah`, so you have to do `[dx,dy] = gradient(blah); [ddx ~] = gradient(dx); [~ ddy] = gradient(dy);`. – Mad Physicist Feb 12 '16 at 17:08
  • @MadPhysicist I am still confused about your Approach for calculating the second derivative. I just looked at this [link](http://de.mathworks.com/matlabcentral/newsreader/view_thread/170336). Dont you think i was doing it the right way. ? – Aakash Sehgal Feb 15 '16 at 09:40
  • @MadPhysicist I just tried out calculating the second Gradient with both the methods and i surprised to see that the results are the same in both cases. But what is different is when i do ` [ddx,ddy] = Gradient([dx,dy])`. Could you please explain the difference ? – Aakash Sehgal Feb 15 '16 at 10:02
  • Yes. I will draft an answer that hopefully addresses all your additional questions. – Mad Physicist Feb 16 '16 at 01:12

1 Answers1

1

To answer your literal question, you should store the variables in a structure array or at least a cell array. If all of your structures have the same fields, you can access all of them by indexing a single array variable, say Structure_element:

for i = 1:numel(Structure_element)
    field = Structure_element(i).value
    % compute gradients of field
end

Now to address the issue of the actual gradient computation. The gradient function computes an approximation for \frac{\partial F}{\partial x}, \frac{\partial F}{\partial y}, where F is your matrix of data. Normally, a MATLAB function is aware of how many output arguments are requested. When you call gradient(gradient(F)), the outer gradient is called on the first output of the inner gradient call. This means that you are currently getting an approximation for \frac{\partial^2 F}{\partial x^2}, \frac{\partial}{\partial y} \frac{\partial F}{\partial x}.

I suspect that you are really trying to get \frac{\partial^2 F}{\partial x^2}, \frac{\partial^2 F}{\partial y^2}. To do this, you have to get both outputs from the inner call to gradient, pass them separately to the outer call, and choose the correct output:

[dx,dy] = gradient(F);
[ddx, ~] = gradient(dx);
[~, ddy] = gradient(dy);

Note the separated calls. The tilde was introduced as a way to ignore function arguments in MATLAB Release 2009b. If you have an older version, just use an actual variable named junk or something like that.

Mad Physicist
  • 107,652
  • 25
  • 181
  • 264