Skip to main content
M
MoeSzyslak12-Amethyst
12-Amethyst
February 5, 2018
Solved

Varying Number of Constraints in a Solve Block (Snell's Law, Mathcad 15)

  • February 5, 2018
  • 5 replies
  • 7261 views

Hi, All --

I'm trying to build a solve block with the following constraints:

 

s1*sin(q1) = s2*sin(q2) = ... = sn*sin(qn) 

 

However, the number of terms (n) is not know a priori. Is there a way to do this in MC15?

 

Thanks.

Matt

Best answer by Werner_E

OK, I couldn't resits 😉

Here's my approach. Looks like the values are close enough to the one you expected and all conditions should be fulfilled. You may even use "Find" instead of "MinErr".

The use of the predefined guess vector for the angles isn't really elegant but the alternative would be to create a suitable guess vector of the correct size before calling the solve block function and add it as additional function argument.

I guess you are only interested in the velocities and not in the various values of the angles in each iteration step. If you should need all the angles, it should be easy enough to collect them in a matrix and return them as well.

Bild1.png

 

BTW, is it possible that the time values you posted are a bit inaccurate 😉

Bild2.png

5 replies

L
LucMeekes23-Emerald IV
23-Emerald IV
February 5, 2018

Is it more difficult than:

LM_20180205_solve.png

 

Success!
Luc

    M
    MoeSzyslak12-AmethystAuthor
    12-Amethyst
    February 5, 2018

    Thanks, Luc -- but that doesn't seem to work -- unless I am doing something else wrong, of course. Please see attached.

    F
    Fred_Kohlhepp23-Emerald I
    23-Emerald I
    February 5, 2018
    wrote:

    Hi, All --

    I'm trying to build a solve block with the following constraints:

     

    s1*sin(q1) = s2*sin(q2) = ... = sn*sin(qn) 

     

    However, the number of terms (n) is not know a priori. Is there a way to do this in MC15?

     

    Thanks.

    Matt

    You haven't supplied enough information.  If you're after the s values (s1, s2, etc) and you know the q's, then Luc has given you an answer.  If you know the s's and are looking for the q's, you can use his answer but's going to take some effort.  If you're after both the s's and q's you need more constraints. . .

    M
    MoeSzyslak12-AmethystAuthor
    12-Amethyst
    February 5, 2018

    Fred -- you're absolutely correct -- I just posted a sheet in response to Luc's solution. I would appreciate it if you could have a look.

     

    Matt

    F
    Fred_Kohlhepp23-Emerald I
    23-Emerald I
    February 5, 2018

    Hmmm

    W
    Werner_E25-Diamond I
    25-Diamond I
    February 5, 2018

    Haven't looked at the last sheet you posted and just tried to make the solve block in your first sheet working.

    Not sure if the following would do the job. I am not proud of this crude solution and hope someone is able to find a more elegant solution.

    Bild.png

    M
    MoeSzyslak12-AmethystAuthor
    12-Amethyst
    February 5, 2018

    Thanks for having a look Werner -- pretty slick way to get the solve block working. Unfortunately, the solution you obtain is quite far from the true solution, shown below. What's odd is that the first set of constraints (the sum) is well-satisfied, but the second (the equality) is not.

    image.png

    M
    MoeSzyslak12-AmethystAuthor
    12-Amethyst
    February 6, 2018

    OK, I was trying to take some shortcuts with this, but perhaps it's best to fully explain the problem. Attached I have produced the iterative solution (brute force) for the first three layers. The solution is, as best I can tell, correct. However, I am looking for a way to formulate it such that I can pass it vectors of layer thicknesses and travel times of arbitrary length -- clear the attached approach is not amenable to that.

     

    Thanks again, all.

    F
    Fred_Kohlhepp23-Emerald I
    23-Emerald I
    February 6, 2018

    But you're changing the problem each time.

    • You solve for the velocity of the first layer:  V1 = sqrt(Z1^2+S^2)/T1  Implying that Z1 is the horizontal distance for a 3 meter thick layer.
    • Then your first solve block says that the sum of z tan(theta) equals the vertical height.  But all of your Z's are 1 meter  (Should be z* cot(theta)   ?)  So now you have the same three meter vertical thickness and two meters horizontal travel.  Not the same problem!
    M
    MoeSzyslak12-AmethystAuthor
    12-Amethyst
    February 6, 2018

    Hi, Fred et al. --

     

    The layer thickness is 1 m, but the horizontal offset of the signal source is 3 m. In general, of course, these values are arbitrary.

     

    I have attached an image of the problem geometry. The problem isn't changing, I am just not going a good job explaining it. I had hoped that I could solve it in one way when I started this thread, but now I am pretty convinced that way won't work. So I have re-cast it a couple of different ways in hopes of better explaining it, but I have served only to sow confusion.

     

    In the attached figure, We have a wave source on the surface and some receivers on a vertical line some horizontal distance away, S. The travel times from source to receivers are measured. In the first layer (V_1), the wave velocity (or slowness) may be calculated directly. That velocity is then used to calculate the wave velocity in layer two (V_2), which requires the knowledge of the refraction angles q_12 and q_22. Thus, the problem is now iterative. Once V_2 is known, the calculation is continued for subsequent layers. In the figure, X_ij corresponds to a property of the jth raypath in the ith layer where X is some quantity (length, angle, etc.).

     

    Does this make sense?

     

    Thanks.

     

    MattCapture.PNG

    W
    Werner_E25-Diamond IAnswer
    25-Diamond I
    February 7, 2018

    OK, I couldn't resits 😉

    Here's my approach. Looks like the values are close enough to the one you expected and all conditions should be fulfilled. You may even use "Find" instead of "MinErr".

    The use of the predefined guess vector for the angles isn't really elegant but the alternative would be to create a suitable guess vector of the correct size before calling the solve block function and add it as additional function argument.

    I guess you are only interested in the velocities and not in the various values of the angles in each iteration step. If you should need all the angles, it should be easy enough to collect them in a matrix and return them as well.

    Bild1.png

     

    BTW, is it possible that the time values you posted are a bit inaccurate 😉

    Bild2.png

      M
      MoeSzyslak12-AmethystAuthor
      12-Amethyst
      February 7, 2018

      This is fantastic, Werner -- very ingenious.

       

      Yes, the time values were imperfect -- just some synthetic numbers to test the algorithm. I'm going to give it a go with some real data.

       

      On a related note, there do exist algorithms that accept not just absolute travel times, but also standard deviations. These are true inversions with regularization, obviously, Those approaches are described in the literature often in a "hand wavy" fashion, making them difficult to implement if you don't want to start at the beginning....

      Terms & ConditionsAccessibility statement