How to loop through all values between a set of values in nested loops and get all combinations of values in python

My goal is to find the minimum of some function "V" on a condition curve "R=0" V and R are both functions of variables X, Y, D, L. I need as outputs all of {X,Y,D,L,V,R} so i can make a design around these parameters (which correspond to physical measurements and values). V and R are too complex to find the minimum analytically so I am trying to solve it numerically. I am using python and would prefer not to change language if I don't need to.

I am trying to loop through all combinations of X, Y, D, & L between some boundaries with small gradations of each (dX, dY, dD, dL) then run each iteration of that through R and V to find the smallest value of V within an acceptable range around R (see below). I am a mechanichal engineering student so apologies if this is a trivial question or poorly worded.

What i currently have is of the simplified form:

def compute_R (X,Y,D,L):
    ...
    return R = compute_R (X,Y,D,L)

def compute_V (X,Y,D,L):
    ...
    return V = compute_V (X,Y,D,L)

def compute_N (X,Y):
    ...
    return N = compute_N (X,Y)

while N > 1: #X loop#
    X += dX
    N = compute_N (X,Y)
    while N > 1: #Y loop#
        Y += dY
        N = compute_N (X,Y)
        while D/2<X and D/2<y: #D loop#
             D += dD
             while L<0.02: #L loop#
                 R = compute_R (X,Y,D,L) 
                 if a<R and b>R:
                     
                     if V < Smallest_V_Yet
                         Smallest_V_Yet=V
                         Best_Solution_Yet= (X,Y,L,D,R,V)

where a and b are values acceptably close to 0; X,Y,L have domains of about [0,0.02], D has a domain of about [0.001,0.02]; dX, dY, dD, dL all have size about 0.0001; all variables are floats, except for Best_Solution_Yet which is a tuple; all of the "..." represent some mathematical relationship, involving intermediate functions and other variables which are not optimization parameters.

The idea is to loop through this system for all combinations of X,Y,D,L in the domains above overwriting Best_Solution_Yet repeatedly until it has iterated through all valid inputs at which point Best_Solution_Yet would be acceptably close to the actual best solution (actual input set which minimizes V and lies exactly on R=0)

This solution, however only only loops through a small number of the solutions that I want to test (it goes to the last while loop then once it gets to L=0.02 it goes to the 2nd to last and so on but i want it to go to the L loop, then go down to the D loop, change D by dD, then go back to L and loop through all valid values of L and propagating in that way).

I have also tried using for loops but I can't get them to accept the conditions that i want them to.

I hope I have been sufficiently clear about my issue. Some final notes are that:

• R, N, and V are themselves composed of other functions which perform functions such as solving the convergence for recursively defined inputs to R (R=1/a*b+c where a=f(b,X,...) and b=f(a,X,...) such that a and b must be simultaneously solved by taking a guess for a=a1, finding b(a1)=b1, a(b1)=a2, b2=b(a2) and so on until abs(1-bn/bn-1)<c and abs(1-an/an-1)<c for some convergence criteria c<<1) and finding intermediate values that are inputs to multiple other inputs to R, N, and/or V.
• I would like everything to be parametric if possible (that is having all of X,Y,D,L,dX,dY,dD,dL,c,R,N,V free to be redefined as necessary)
• This does not need to be extremely computationally efficient, however, any advice for improving the speed would be appreciated (other than not using python)

Edit: incorporating formatting feedback