Skip to content

Commit

Permalink
Commitin
Browse files Browse the repository at this point in the history
  • Loading branch information
dafeda committed Aug 6, 2024
1 parent 285e686 commit 8024730
Show file tree
Hide file tree
Showing 3 changed files with 40 additions and 38 deletions.
32 changes: 32 additions & 0 deletions test-data/heat_equation/definition.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
from collections import namedtuple

import numpy as np

# Number of grid-cells in x and y direction
nx = 10

# time steps
k_start = 0
k_end = 500

# Define initial condition, i.e., the initial temperature distribution.
# How you define initial conditions will effect the spread of results,
# i.e., how similar different realisations are.
u_init = np.zeros((k_end, nx, nx))
u_init[:, 5:7, 5:7] = 100

# Resolution in the x-direction (nothing to worry about really)
dx = 1

Coordinate = namedtuple("Coordinate", ["x", "y"])

obs_coordinates = [
Coordinate(5, 3),
Coordinate(3, 5),
Coordinate(5, 7),
Coordinate(7, 5),
Coordinate(2, 2),
Coordinate(7, 2),
]

obs_times = np.linspace(10, k_end, 8, endpoint=False, dtype=int)
15 changes: 1 addition & 14 deletions test-data/heat_equation/generate_files.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,14 +2,14 @@
Contains code that was used to generate files expected by ert.
"""

from collections import namedtuple
from textwrap import dedent
from typing import Callable, List

import numpy as np
import numpy.typing as npt
import pandas as pd
import xtgeo
from definition import Coordinate, obs_coordinates, obs_times
from heat_equation import heat_equation, sample_prior_conductivity

rng = np.random.default_rng(11223344)
Expand All @@ -23,9 +23,6 @@ def create_egrid_file():
grid.to_file("CASE.EGRID", "egrid")


Coordinate = namedtuple("Coordinate", ["x", "y"])


def make_observations(
coordinates: List[Coordinate],
times: npt.NDArray[np.int_],
Expand Down Expand Up @@ -101,16 +98,6 @@ def make_observations(
dt = dx**2 / (4 * max(np.max(cond_truth), np.max(cond_truth)))

u_t = heat_equation(u_init, cond_truth, dx, dt, k_start, k_end, rng=rng)
obs_coordinates = [
Coordinate(5, 3),
Coordinate(3, 5),
Coordinate(5, 7),
Coordinate(7, 5),
Coordinate(2, 2),
Coordinate(7, 2),
]

obs_times = np.linspace(10, k_end, 8, endpoint=False, dtype=int)

d = make_observations(
obs_coordinates, obs_times, u_t, lambda value: abs(0.05 * value)
Expand Down
31 changes: 7 additions & 24 deletions test-data/heat_equation/heat_equation.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,11 +3,12 @@

from typing import Optional

import geostat
import numpy as np
import numpy.typing as npt
from definition import dx, k_end, k_start, nx, obs_coordinates, obs_times, u_init

rng = np.random.default_rng(1234)
import geostat
import numpy.typing as npt


def heat_equation(
Expand Down Expand Up @@ -57,22 +58,6 @@ def sample_prior_conductivity(ensemble_size, nx, rng):


if __name__ == "__main__":
# Number of grid-cells in x and y direction
nx = 10

# time steps
k_start = 0
k_end = 500

# Define initial condition, i.e., the initial temperature distribution.
# How you define initial conditions will effect the spread of results,
# i.e., how similar different realisations are.
u_init = np.zeros((k_end, nx, nx))
u_init[:, 5:7, 5:7] = 100

# Resolution in the x-direction (nothing to worry about really)
dx = 1

cond = sample_prior_conductivity(ensemble_size=1, nx=nx, rng=rng).reshape(nx, nx)

# Write the array to a GRDECL formatted file
Expand All @@ -95,9 +80,7 @@ def sample_prior_conductivity(ensemble_size, nx, rng):

response = heat_equation(u_init, cond, dx, dt, k_start, k_end, rng)

# Pick the responses where you have observations
index = [0, 2, 5]
with open("gen_data_0.out", "w", encoding="utf-8") as f:
for i in index:
for value in response[0][i]:
f.write(f"{value:.6f}\n")
index = sorted((obs.x, obs.y) for obs in obs_coordinates)
for time_step in obs_times:
with open(f"gen_data_{time_step}.out", "w", encoding="utf-8") as f:
f.write("\\n".join(str(response[time_step][i]) for i in index))

0 comments on commit 8024730

Please sign in to comment.