feature: combine multiple assimilation methods

data_assimilation/ARTSS.py

@@ -1,11 +1,61 @@
 #!/usr/bin/env python3
 import os
 import xml.etree.ElementTree as ET
-from typing import List, Dict, Union, Set
+from typing import List, Dict, Union, Set, TextIO
 
 import numpy as np
 
 from obstacle import Obstacle, PATCHES, STATE
+from utility import log
+
+
+def start_new_instance(output_file: str, directory: str, artss_exe_path: str,
+                       artss_exe: str = 'artss_data_assimilation_serial'):
+    cwd = os.getcwd()
+    os.chdir(directory)
+    exe_command = f'mpirun --np 2 {os.path.join(artss_exe_path, artss_exe)} {output_file}'
+    print(os.getcwd(), exe_command)
+    os.system(exe_command)
+    os.chdir(cwd)
+
+
+def change_xml(change: Dict[str, float], input_file: str, output_file: str, artss_root_path: str, artss_data_path: str, file_debug: TextIO):
+    out_file = os.path.join(artss_data_path, output_file)
+    log(f'out_file: {out_file}', file_debug)
+
+    command = ''
+    for key in change:
+        command += f' --{key} {change[key]}'
+    command = f'{os.path.join(artss_root_path, "tools", "change_xml.sh")} -i {input_file} -o {out_file} --loglevel off' + command
+    log(f'{command}', file_debug)
+    os.system(command)
+
+
+def create_start_xml(change: dict, input_file: str, output_file: str, t_revert: float, file: str, t_end: float,
+                     directory: str, artss_path: str, file_debug: TextIO, dir_name: str = '') -> [str, str]:
+    f = os.path.join('..', file)
+    f = f.replace("/", "\/")
+
+    command = ''
+    name = ''
+    for key in change:
+        command += f' --{key} {change[key]}'
+        name += f'{key}_'
+    name = name[:-1]
+
+    if not dir_name == '':
+        name = dir_name
+    output_dir = os.path.join(directory, '.vis', name)
+    if not os.path.exists(output_dir):
+        os.makedirs(output_dir)
+
+    out_file = os.path.join(output_dir, output_file)
+    log(f'out_file: {out_file}', file_debug)
+    log(f'out_dir: {output_dir}', file_debug)
+    command = f'{os.path.join(artss_path, "tools", "change_xml.sh")} -i {input_file} -o {out_file} --da {t_revert} "{f}" 7779 --tend {t_end} --loglevel off' + command
+    log(f'{command}', file_debug)
+    os.system(command)
+    return out_file, output_dir
 
 
 class XML:
@@ -52,7 +102,7 @@ def get_output_time_step(self) -> float:
         return 1
 
     def get_temperature_source(self) -> dict:
-        if self.temperature_source is None:
+        if not self.temperature_source:
             # check if temperature source is even there ? or crash
             root = self.xml_tree.getroot()
             source_tree = root.find('solver').find('temperature').find('source')
@@ -89,8 +139,7 @@ def __init__(self, domain_param: Dict[str, float], enable_computational_domain:
         for o in obstacles:
             self.obstacles[o.name] = o
         self.domain_boundary_cells: Dict[str, Set[int]] = {}
-        self.domain_inner_cells: Dict[str, Set[int]] = {}
-        self.obstacle_cells: Dict[str, Dict[str, Set[int]]] = {}
+        self.domain_inner_cells: Set[int] = set()
         self.computational_domain(domain_param, enable_computational_domain)
         self.calculate_domain()
         self.calculate_indices()
@@ -166,12 +215,12 @@ def calculate_obstacle_cells(self, o: Obstacle):
         index_y2: int = self.get_matching_obstacle_index(o.geometry['oy2'], 'y')
         index_z1: int = self.get_matching_obstacle_index(o.geometry['oz1'], 'z') + 1
         index_z2: int = self.get_matching_obstacle_index(o.geometry['oz2'], 'z')
+        o.index = {'x1': index_x1, 'x2': index_x2, 'y1': index_y1, 'y2': index_y2, 'z1': index_z1, 'z2': index_z2}
         for i in range(index_x1 + 1, index_x2):
             for j in range(index_y1 + 1, index_y2):
                 for k in range(index_z1 + 1, index_z2):
                     inner_cells.append(self.calculate_index(i, j, k))
-        self.obstacle_cells[o.name]: Dict[str, Set[int]] = {}
-        self.obstacle_cells[o.name]['inner cells'] = set(inner_cells)
+        self.obstacles[o.name].cells['inner cells'] = set(inner_cells)
 
         # indices of obstacle boundary cells
         front = []
@@ -180,26 +229,26 @@ def calculate_obstacle_cells(self, o: Obstacle):
             for j in range(index_y1, index_y2 + 1):
                 front.append(self.calculate_index(i, j, index_z1))
                 back.append(self.calculate_index(i, j, index_z2))
-        self.obstacle_cells[o.name]['front'] = set(front)
-        self.obstacle_cells[o.name]['back'] = set(back)
+        self.obstacles[o.name].cells['front'] = set(front)
+        self.obstacles[o.name].cells['back'] = set(back)
 
         bottom = []
         top = []
         for i in range(index_x1, index_x2 + 1):
             for k in range(index_z1, index_z2 + 1):
                 bottom.append(self.calculate_index(i, index_y1, k))
                 top.append(self.calculate_index(i, index_y2, k))
-        self.obstacle_cells[o.name]['bottom'] = set(bottom)
-        self.obstacle_cells[o.name]['top'] = set(top)
+        self.obstacles[o.name].cells['bottom'] = set(bottom)
+        self.obstacles[o.name].cells['top'] = set(top)
 
         left = []
         right = []
         for j in range(index_y1, index_y2 + 1):
             for k in range(index_z1, index_z2 + 1):
                 left.append(self.calculate_index(index_x1, j, k))
                 right.append(self.calculate_index(index_x2, j, k))
-        self.obstacle_cells[o.name]['left'] = set(left)
-        self.obstacle_cells[o.name]['right'] = set(right)
+        self.obstacles[o.name].cells['left'] = set(left)
+        self.obstacles[o.name].cells['right'] = set(right)
 
     def calculate_domain_boundaries(self):
         front = []
@@ -236,9 +285,9 @@ def calculate_domain_inner(self):
                 for k in range(self.domain_param['start z index CD'], self.domain_param['end z index CD']):
                     inner.append(self.calculate_index(i, j, k))
         inner = set(inner)
-        for oc in self.obstacle_cells:
-            for type_of_cell in self.obstacle_cells[oc]:
-                inner = inner - self.obstacle_cells[oc][type_of_cell]
+        for obstacle in self.obstacles.values():
+            for type_of_cell in obstacle.cells:
+                inner = inner - obstacle.cells[type_of_cell]
         self.domain_inner_cells = inner
 
     def calculate_index(self, i, j, k):
@@ -250,14 +299,14 @@ def print_info(self):
             f"Domain size inner cells: {self.domain_param['nx']} {self.domain_param['ny']} {self.domain_param['nz']}\n"
             f"step size (x|y|z): ({self.domain_param['dx']}|{self.domain_param['dy']}|{self.domain_param['dz']})")
         for o in self.obstacles.values():
-            print(f"-- Obstacle {o.name}\n"
+            print(f"-- Obstacle {o.name} ({o.state})\n"
                   f"   size of slices (Front|Back Bottom|Top Left|Right): "
-                  f"{len(self.obstacle_cells[o.name]['front'])}|"
-                  f"{len(self.obstacle_cells[o.name]['back'])} "
-                  f"{len(self.obstacle_cells[o.name]['bottom'])}|"
-                  f"{len(self.obstacle_cells[o.name]['top'])} "
-                  f"{len(self.obstacle_cells[o.name]['left'])}|"
-                  f"{len(self.obstacle_cells[o.name]['right'])}\n"
+                  f"{len(self.obstacles[o.name].cells['front'])}|"
+                  f"{len(self.obstacles[o.name].cells['back'])} "
+                  f"{len(self.obstacles[o.name].cells['bottom'])}|"
+                  f"{len(self.obstacles[o.name].cells['top'])} "
+                  f"{len(self.obstacles[o.name].cells['left'])}|"
+                  f"{len(self.obstacles[o.name].cells['right'])}\n"
                   f"   coords (x y z): "
                   f"({o.geometry['ox1']}|{o.geometry['ox2']}) "
                   f"({o.geometry['oy1']}|{o.geometry['oy2']}) "
@@ -279,23 +328,27 @@ def print_debug(self):
         print(f"index X PD: ({self.domain_param['start x index PD']}|{self.domain_param['end x index PD']}) "
               f"index Y PD: ({self.domain_param['start y index PD']}|{self.domain_param['end y index PD']}) "
               f"index Z PD: ({self.domain_param['start z index PD']}|{self.domain_param['end z index PD']})")
+        print(f"inner cells index: ({min(self.domain_inner_cells)}|{max(self.domain_inner_cells)})")
+        for patch in PATCHES:
+            print(f"boundary cells index {patch}: ({min(self.domain_boundary_cells[patch])}|{max(self.domain_boundary_cells[patch])})")
+        print()
         for o in self.obstacles.values():
             print(f"-- Obstacle {o.name}\n"
                   f"   size of slices (Front|Back Bottom|Top Left|Right): "
-                  f"{len(self.obstacle_cells[o.name]['front'])}|"
-                  f"{len(self.obstacle_cells[o.name]['back'])} "
-                  f"{len(self.obstacle_cells[o.name]['bottom'])}|"
-                  f"{len(self.obstacle_cells[o.name]['top'])} "
-                  f"{len(self.obstacle_cells[o.name]['left'])}|"
-                  f"{len(self.obstacle_cells[o.name]['right'])}\n"
+                  f"{len(self.obstacles[o.name].cells['front'])}|"
+                  f"{len(self.obstacles[o.name].cells['back'])} "
+                  f"{len(self.obstacles[o.name].cells['bottom'])}|"
+                  f"{len(self.obstacles[o.name].cells['top'])} "
+                  f"{len(self.obstacles[o.name].cells['left'])}|"
+                  f"{len(self.obstacles[o.name].cells['right'])}\n"
                   f"   coords (x y z): "
                   f"({o.geometry['ox1']}|{o.geometry['ox2']}) "
                   f"({o.geometry['oy1']}|{o.geometry['oy2']}) "
                   f"({o.geometry['oz1']}|{o.geometry['oz2']})")
             str_obstacle_cells = ''
             for patch in PATCHES:
-                start = min(self.obstacle_cells[o.name][patch])
-                end = max(self.obstacle_cells[o.name][patch])
+                start = min(self.obstacles[o.name].cells[patch])
+                end = max(self.obstacles[o.name].cells[patch])
                 str_obstacle_cells += f"   ({patch}: {start}|{end})\n"
                 k = self.get_k(start)
                 j = self.get_j(start, k=k)
@@ -305,7 +358,7 @@ def print_debug(self):
                 j = self.get_j(end, k=k)
                 i = self.get_i(end, k=k, j=j)
                 str_obstacle_cells += f"   {patch} end: {i}|{j}|{k}\n"
-            print(str_obstacle_cells)
+            print(str_obstacle_cells[:-2])
 
     def get_ijk_from_xyz(self, coord_x, coord_y, coord_z):
         return self.get_i_from_x(coord_x), self.get_j_from_y(coord_y), self.get_k_from_z(coord_z)
@@ -336,11 +389,11 @@ def get_type(self, index):
         for p in PATCHES.keys():
             if index in self.domain_boundary_cells[p]:
                 matches.append(f'domain boundary cell ({p})')
-        for oc in self.obstacle_cells:
+        for oc in self.obstacles:
             for p in PATCHES.keys():
-                if index in self.obstacle_cells[oc][p]:
+                if index in self.obstacles[oc].cells[p]:
                     matches.append(f'obstacle boundary cell ({oc}/{p})')
-            if index in self.obstacle_cells['inner cells']:
+            if index in self.obstacles[oc].cells['inner cells']:
                 matches.append(f'obstacle inner cell ({oc})')
         if len(matches) == 0:
             return "cell type unknown"
@@ -356,17 +409,22 @@ def change_obstacle(self, obstacle: Obstacle):
         self.obstacles[obstacle.name] = obstacle
         self.calculate_obstacle_cells(obstacle)
 
-    def remove_obstacle(self, name: str):
-        self.obstacles[name].state = 'deleted'
+    def remove_obstacle(self, name: str) -> Obstacle:
+        obstacle = self.obstacles.pop(name)
+        obstacle.state = 'deleted'
+        return obstacle
 
     def set_value_of_obstacle_cells(self, value: float, field: np.ndarray, obstacle_name: str):
         cells: Set[float] = set()
-        for val in self.obstacle_cells[obstacle_name].values():
+        for val in self.obstacles[obstacle_name].cells.values():
             cells.update(val)
 
         for cell in cells:
             field[cell] = value
 
     def set_value_of_obstacle_patch(self, value: float, field: np.ndarray, obstacle_name: str, patch: str):
-        for cell in self.obstacle_cells[obstacle_name][patch]:
+        for cell in self.obstacles[obstacle_name].cells[patch]:
             field[cell] = value
+
+    def get_obstacles(self):
+        return list(self.obstacles.values())

data_assimilation/TCP_client.py

@@ -1,11 +1,14 @@
+import datetime
 import socket
+import struct
 
 
 class TCPClient:
     def __init__(self):
         # Create a TCP/IP socket
         self.socket = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
-
+        timeval = struct.pack('ll', 20, 0)
+        self.socket.setsockopt(socket.SOL_SOCKET, socket.SO_SNDTIMEO, timeval)
         # Connect the socket to the port where the server is listening
         self.server_address = ('localhost', 7777)
 
@@ -23,7 +26,8 @@ def send_message(self, message: bin):
         # Look for the response
         expected_response = "message was received: Rollback done"
         response = ''
-
+        self.socket.settimeout(30)
+        print('time', datetime.datetime.now())
         while len(response) < len(expected_response):
             try:
                 reply = self.socket.recv(16)