LoopStructural/tests/unit/interpolator/test_discrete_supports.py at 7c133732c8a8eb50c77d3605bbb8124c32615915 · Loop3D/LoopStructural · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
from LoopStructural.interpolators import StructuredGrid
import numpy as np
import pytest

## structured grid tests
def test_create_support(support):
    """
    support is a fixture that returns a support object created
    with the default constructor. Ensure that it is not none and
    make sure the origin and maximum are correct.
    """
    assert support is not None
    assert np.sum(support.origin - np.zeros(3)) == 0
    assert np.sum(support.maximum - np.ones(3) * 10) == 0


def test_create_support_origin_nsteps(support_class):
    grid = support_class(
        origin=np.zeros(3),
        nsteps=np.array([10, 10, 10]),
        step_vector=np.array([0.1, 0.1, 0.1]),
    )
    assert np.sum(grid.step_vector - np.array([0.1, 0.1, 0.1])) == 0
    assert np.sum(grid.maximum - np.ones(3)) == 0


def test_inside(support):
    assert np.all(support.inside(support.barycentre))


def test_evaluate_value(support):
    # print()
    # print(support.evaluate_value(support.barycentre, support.nodes[:, 0]))
    # print(support.barycentre[:, 0])
    assert (
        np.sum(
            support.barycentre[:, 0]
            - support.evaluate_value(support.barycentre, support.nodes[:, 0])
        )
        == 0
    )

@pytest.mark.parametrize('steps',[10,20,100])
def test_evaluate_gradient(support_class,steps):
    support = support_class(nsteps=[steps]*3)
    # test by setting the scalar field to the y coordinate
    vector = support.evaluate_gradient(support.barycentre, support.nodes[:, 1])
    assert np.sum(vector - np.array([0, 1, 0])) == 0

    # # same test but for a bigger grid, making sure scaling for cell is ok
    # support = support_class(step_vector=np.array([100, 100, 100]),)
    # vector = support.evaluate_gradient(support.barycentre, support.nodes[:, 1])
    # assert np.sum(vector - np.array([0, 1, 0])) == 0


def test_outside_box(support):
    # test by setting the scalar field to the y coordinate
    inside = support.inside(support.barycentre + 5)
    assert np.all(~inside == np.any((support.barycentre + 5) > support.maximum, axis=1))
    inside = support.inside(support.barycentre - 5)
    assert np.all(~inside == np.any((support.barycentre - 5) < support.origin, axis=1))

    cell_indexes, inside = support.position_to_cell_index(support.barycentre - 5)
    assert np.all(cell_indexes[inside, 0] < support.nsteps_cells[0])
    assert np.all(cell_indexes[inside, 1] < support.nsteps_cells[1])
    assert np.all(cell_indexes[inside, 2] < support.nsteps_cells[2])
    corners = support.cell_corner_indexes(cell_indexes)
    assert np.all(corners[inside, 0] < support.nsteps[0])
    assert np.all(corners[inside, 1] < support.nsteps[1])
    assert np.all(corners[inside, 2] < support.nsteps[2])
    globalidx = support.global_node_indices(corners)
    # print(globalidx[inside],grid.n_nodes,inside)
    assert np.all(globalidx[inside] < support.n_nodes)
    inside = support.inside(support.barycentre - 5)
    # inside, support.position_to_cell_corne    rs(support.barycentre - 5)
    vector = support.evaluate_gradient(support.barycentre - 5, support.nodes[:, 1])
    assert np.sum(np.mean(vector[inside, :], axis=0) - np.array([0, 1, 0])) == 0
    vector = support.evaluate_gradient(support.nodes, support.nodes[:, 1])


@pytest.mark.parametrize("seed", range(10))
def test_evaluate_gradient2(support_class, seed):
    rng = np.random.default_rng(seed)
    step = rng.uniform(0, 100)
    grid = support_class(step_vector=np.array([step, step, step]))

    # define random vector
    n = rng.random(3)
    n /= np.linalg.norm(n)
    distance = n[0] * grid.nodes[:, 0] + n[1] * grid.nodes[:, 1] + n[2] * grid.nodes[:, 2]
    vector = grid.evaluate_gradient(rng.uniform(1, 8, size=(100, 3)), distance)
    assert np.all(np.isclose(np.sum(vector - n[None, :], axis=1), 0, atol=1e-3, rtol=1e-3))


def test_get_element(support):
    point = support.barycentre[[0], :]
    # point[0, 0] += 0.1
    vertices, dof, idc, inside = support.get_element_for_location(point)
    # vertices = vertices.reshape(-1, 3)
    bary = np.mean(vertices, axis=1)
    assert np.isclose(np.sum(point - bary), 0)


def test_global_to_local_coordinates():
    grid = StructuredGrid()
    point = np.array([[1.2, 1.5, 1.7]])
    local_coords = grid.position_to_local_coordinates(point)
    assert np.isclose(local_coords[0, 0], 0.2)
    assert np.isclose(local_coords[0, 1], 0.5)
    assert np.isclose(local_coords[0, 2], 0.7)


def test_get_element_outside(support):
    point = np.array([support.origin - np.ones(3)])
    idc, inside = support.position_to_cell_corners(point)
    assert not inside[0]


def test_node_index_to_position(support):
    assert (
        np.sum(support.node_indexes_to_position(np.array([[0, 0, 0]])) - np.array([0, 0, 0])) == 0
    )
    for i in range(10):
        for j in range(10):
            for k in range(10):
                assert (
                    np.sum(
                        support.node_indexes_to_position(np.array([[i, j, k]]))
                        - np.array([i, j, k]) * support.step_vector
                    )
                    == 0
                )
    assert np.sum(support.node_indexes_to_position(np.array([0, 0, 0])) - np.array([0, 0, 0])) == 0


def test_global_index_to_cell_index(support):
    assert np.sum(support.global_index_to_cell_index(np.array([0])) - np.array([0, 0, 0])) == 0


def test_global_index(support):
    indexes = np.array(
        np.meshgrid(
            np.arange(0, support.nsteps[0]),
            np.arange(0, support.nsteps[1]),
            np.arange(0, support.nsteps[2]),
        )
    ).reshape(-1, 3)
    global_node_index = support.global_node_indices(indexes)
    assert np.all(global_node_index >= 0)
    assert np.all(global_node_index < support.n_nodes)

    indexes = np.array(
        np.meshgrid(
            np.arange(0, 3),
            np.arange(0, 1),
            np.arange(0, 1),
        )
    ).reshape(-1, 3)
    global_node_index = support.global_node_indices(indexes)
    assert np.all(global_node_index >= 0)
    assert np.all(global_node_index < support.n_nodes)