Re: how to treat different kind of errors when using numeric solvers?

Sorry, i cut of part of the function.


def detect_cable(data: np.ndarray):
    b_measured = np.average(data, axis=1)
    sensor_coordinates = np.array([[.265, 0, .382], [0, .712, .764], [0, .712, 0], [.752, .712, .382]])
    sensor_num = b_measured.shape[0]

    # print(b_measured)
    guess_b_earth = np.array(
        [1.7921, 4.6676, -1.013])  # https://ngdc.noaa.gov/geomag/calculators/magcalc.shtml#igrfwmm
    # Latitude:    54° 08' 47" N
    # Longitude:   8° 48' 50" E
    guess_n = np.array([0, 0, 1])
    guess_a = np.array([.2, 0, 0])
    guess_r = np.array([0, 2.7, 0])
    guess_i = 1000
    p0 = np.hstack([guess_b_earth, guess_n, guess_r, guess_a, guess_i])

    def error_func_geometry_vs_measurement(x, b_measured):
        b_earth = x[0:3]
        n = x[3:6]
        r = x[6:9]
        a = x[9:12]
        i = x[12]

        k = i * MU_0 / np.pi
        for sensor_cnt in range(0, sensor_num):
            r_relativ = 0
            if sensor_cnt:  # only do this for sensor 1,2 and 3, not 0.
                # that means that the algorithm will find the projection of 0's point on to the cable on its own.
                # The other r_relativ vectors depend on r
                r_relativ = r - n * np.matmul((r - sensor_coordinates[sensor_cnt]), n) / np.linalg.norm(n)
            else: # for sensor 0 we take r as the projection point on to the cable
                r_relativ = r - sensor_coordinates[sensor_cnt]
            r_magnitude = np.linalg.norm(r_relativ, keepdims=True)
#            print("r: {}".format(r_relativ))
            error_0 = -k * np.cross(a - 2 * (r_relativ.T * np.matmul(r_relativ, a)).T, n) / r_magnitude ** 4 + b_earth - b_measured
            error_0 = np.linalg.norm(error_0, axis=1)
            error_0 = np.dot(error_0, error_0)
#        error_1 = np.dot(n, a) ** 2
#        error_2 = np.linalg.norm(np.matmul(r, n)) ** 2
        error_3 = error_0 # + error_1 + error_2

#        print("error_3 {}, error_0 {}, error_1 {}, error_2 {}".format(error_3, error_0, error_1, error_2))
        return error_3

    res = minimize(error_func_geometry_vs_measurement, p0, args=b_measured)
    x = res.x
    b_earth = x[0:3]
    n = x[3:6]
    r = x[6:9]
    a = x[9:12]
    i = x[12]

    print("B_earth: {}, norm: {}".format(b_earth, np.linalg.norm(b_earth)))
    print("r: {}, Depth: {}, distance {}".format(r, r[1]-sensor_coordinates[1,1], np.linalg.norm(r)))
    print("a {}, norm: {}, dot(a,n) {}".format(a, np.linalg.norm(a), np.dot(a,n)))
    print("n {}, norm: {}, dot(r,n) {}".format(n, np.linalg.norm(a), np.dot(r,n)))
    print("current I: {}".format(i))


________________________________
Von: Schuldei, Andreas
Gesendet: Donnerstag, 11. November 2021 14:49:07
An: [email protected]
Betreff: how to treat different kind of errors when using numeric solvers?


I would like some help from experienced users of numeric solvers (like optimize.minimize()).


I have a system of vector equations that I want to solve numerically.

[cid:6d69c188-653b-4b2b-9f36-69bed2a70a15]
How should I combine the errors, so the solver can use the feedback best? The resulting errors of the first two equations are about 10 orders of magnitude smaller then the third.
What is a "natural" mathematical way to combine such errors?

this is my code:


def detect_cable(data: np.ndarray):
    b_measured = np.average(data, axis=1)
    sensor_coordinates = np.array([[.265, 0, .382], [0, .712, .764], [0, .712, 0], [.752, .712, .382]])
    sensor_num = b_measured.shape[0]

    # print(b_measured)
    guess_b_earth = np.array(
        [1.7921, 4.6676, -1.013])  # https://ngdc.noaa.gov/geomag/calculators/magcalc.shtml#igrfwmm
    # Latitude:    54° 08' 47" N
    # Longitude:   8° 48' 50" E
    guess_n = np.array([0, 0, 1])
    guess_a = np.array([.2, 0, 0])
    guess_r = np.array([0, 2.7, 0])
    guess_i = 1000
    p0 = np.hstack([guess_b_earth, guess_n, guess_r, guess_a, guess_i])

    def error_func_geometry_vs_measurement(x, b_measured):
        b_earth = x[0:3]
        n = x[3:6]
        r = x[6:9]
        a = x[9:12]
        i = x[12]

        k = i * MU_0 / np.pi
        for sensor_cnt in range(0, sensor_num):
            r_relativ = 0
            if sensor_cnt:  # only do this for sensor 1,2 and 3, not 0.
                # that means that the algorithm will find the projection of 0's point on to the cable on its own.
                # The other r_relativ vectors depend on r
                r_relativ = r - n * np.matmul((r - sensor_coordinates[sensor_cnt]), n) / np.linalg.norm(n)
            else: # for sensor 0 we take r as the projection point on to the cable
                r_relativ = r - sensor_coordinates[sensor_cnt]
            r_magnitude = np.linalg.norm(r_relativ, keepdims=True)
#            print("r: {}".format(r_relativ))
            error_0 = -k * np.cross(a - 2 * (r_relativ.T * np.matmul(r_relativ, a)).T, n) / r_magnitude ** 4 + b_earth - b_measured
            error_0 = np.linalg.norm(error_0, axis=1)
            error_0 = np.dot(error_0, error_0)
#        error_1 = np.dot(n, a) ** 2
#        error_2 = np.linalg.norm(np.matmul(r, n)) ** 2
        error_3 = error_0 # + error_1 + error_2

#        print("error_3 {}, error_0 {}, error_1 {}, error_2 {}".format(error_3, error_0, error_1, error_2))
        return error_3

_______________________________________________
SciPy-User mailing list -- [email protected]
To unsubscribe send an email to [email protected]
https://mail.python.org/mailman3/lists/scipy-user.python.org/
Member address: [email protected]