### ===================================================================================================================
### Helper functions for geo-output handling
### ===================================================================================================================
# Copyright ©VIIA 2024
### ===================================================================================================================
### 1. Import modules
### ===================================================================================================================
# General imports
from __future__ import annotations
from typing import TYPE_CHECKING, List
from collections import namedtuple
from pathlib import Path
import json
# References for functions and classes in the rhdhv_fem package
from rhdhv_fem.shape_geometries import Polyline
from rhdhv_fem.shapes import Lines, LineMass, Fstrip
from rhdhv_fem.fem_math import fem_greater, fem_smaller
# References for functions and classes in the viiaPackage
if TYPE_CHECKING:
from viiapackage.viiaStatus import ViiaProject
### ===================================================================================================================
### 2. Helper functions for Geo Output NLTH
### ===================================================================================================================
[docs]def _calculate_forces_per_timestep_fstrip(dat_file, fstrip_itf, time_step_numbers, time_step, result_dict):
""" Function for creating forces for a mapped fstrip per time step."""
Force_data = namedtuple('Force_data', ['base_shear_x', 'base_shear_y', 'base_shear_z', 'timestep'])
# Find all elements for the fstrips by element set
fstrip_eqv_elems = dat_file.elementsets[fstrip_itf]['elements']
# Base shear forces
base_shear_x = 0
base_shear_y = 0
base_shear_z = 0
# Calculating the total forces in x, y, z per strip
for element in fstrip_eqv_elems:
element_intpnt_nr = len((result_dict[str(element)]))
element_area = dat_file.elements[element]['area']
# Get the result for the base shear in all three directions across all integration points
base_shear_x_force = (sum(
[i['STSx'][time_step] for i in result_dict[str(element)].values()]) / element_intpnt_nr) * element_area
base_shear_y_force = (sum(
[i['STSy'][time_step] for i in result_dict[str(element)].values()]) / element_intpnt_nr) * element_area
base_shear_z_force = (sum(
[i['STNz'][time_step] for i in result_dict[str(element)].values()]) / element_intpnt_nr) * element_area
base_shear_x += base_shear_x_force
base_shear_y += base_shear_y_force
base_shear_z += base_shear_z_force
# The force data containing the forces for every time step number is returned. Here the time steps would
# start from 1
return Force_data(
base_shear_x=base_shear_x, base_shear_y=base_shear_y, base_shear_z=base_shear_z,
timestep=time_step_numbers[time_step - (time_step_numbers[0] - 1)] - (time_step_numbers[0] - 1)), \
fstrip_eqv_elems
[docs]def _calculate_forces_per_timestep_piles(pile, time_step_numbers, time_step, result_dict):
""" Function for creating forces for piles per time step."""
Force_data = namedtuple('Force_data', ['base_shear_x', 'base_shear_y', 'base_shear_z', 'timestep'])
# Saving the forces per pile
if pile.pile_node in list(result_dict.keys()):
base_shear_x = result_dict[pile.pile_node]['FNTX'][time_step]
base_shear_y = result_dict[pile.pile_node]['FNTY'][time_step]
base_shear_z = result_dict[pile.pile_node]['FNTZ'][time_step]
else:
raise ValueError(
f"ERROR: The results for the pile node {pile.pile_node} are not available. Please check your output.")
# The force data containing the forces for every time step number is returned. Here the time steps would
# start from 1
return Force_data(
base_shear_x=base_shear_x, base_shear_y=base_shear_y, base_shear_z=base_shear_z,
timestep=time_step_numbers[time_step - (time_step_numbers[0] - 1)] - (time_step_numbers[0] - 1))
[docs]def _get_fstrip_width(p1: List[float], p2: List[float], p3: List[float], p4: List[float]) -> float:
"""
Input:
- p1 (list of 3 float): The coordinate of corner point that defines the quadrilateral.
- p2 (list of 3 float): The coordinate of corner point that defines the quadrilateral.
- p3 (list of 3 float): The coordinate of corner point that defines the quadrilateral.
- p4 (list of 3 float): The coordinate of corner point that defines the quadrilateral.
Output:
- Returns the width of the fstrips, in [m].
"""
# Group randomly until it finds the parallel lines
import math
candi = [p2, p3, p4]
for pt in candi:
# First group
grp_1 = [p1, pt]
grp_2 = [x for x in candi if x is not pt]
# If the slope is inf
if (grp_1[0][0] - grp_1[1][0]) == 0 and (grp_2[0][0] - grp_2[1][0]) == 0:
return math.fabs(grp_1[0][0] - grp_2[0][0])
elif (grp_1[0][0] - grp_1[1][0]) == 0 and (grp_2[0][0] - grp_2[1][0]) != 0:
continue
elif (grp_1[0][0] - grp_1[1][0]) != 0 and (grp_2[0][0] - grp_2[1][0]) == 0:
continue
else:
slope_1 = (grp_1[0][1] - grp_1[1][1]) / (grp_1[0][0] - grp_1[1][0])
slope_2 = (grp_2[0][1] - grp_2[1][1]) / (grp_2[0][0] - grp_2[1][0])
if round(slope_1, 2) == round(slope_2, 2):
break
else:
raise ValueError("ERROR: Parallel lines cannot be found, width cannot be calculated.")
m = slope_1
x1, y1 = grp_1[0][0], grp_1[0][1]
x2, y2 = grp_2[0][0], grp_2[0][1]
d = math.fabs(y1 - y2 - m * x1 + m * x2) / math.sqrt(m ** 2 + 1)
return d
def _get_corner_nodes(shape: Fstrip):
# All shape lines
shape_lines = shape.contour.get_lines()
conner_nodes = []
for node in shape.contour.get_nodes():
node_conn_lines = node.get_lines()
# Get the intersection line
node_lines = [line for line in node_conn_lines if line in shape_lines]
if len(node_lines) != 2:
raise ValueError(f"ERROR: The {node} has more than two lines connected at the contour, please check it.")
# If the left/right lines around the nodes have a different orientation, then the node is in corner
dir_1 = node_lines[0].get_direction().vector
dir_2 = node_lines[1].get_direction().vector
from rhdhv_fem.fem_math import fem_compare_coordinates, fem_flip_vector
if not (fem_compare_coordinates(dir_1, dir_2) or fem_compare_coordinates(dir_1, fem_flip_vector(dir_2))):
conner_nodes.append(node.coordinates)
if len(conner_nodes) != 4:
shape.project.write_log(f"WARNING: There are more than 4 nodes for {shape}, please check, procedure continues.")
return conner_nodes
[docs]def viia_append_envelope_requests_geo_abaqus(
project: ViiaProject, json_path: Path, pile_displacement_requets: bool = False):
"""
Appends json file that is an input for abaqus envelope script with the data about the section forces of the
foundation strips.
Input:
- project (obj): VIIA project object containing collections of fem objects and project variables.
- json_path (Path): Path object of the json to be appended location.
"""
from rhdhv_fem.abaqus_utils.helper_functions.adjust_detailing import _adjust_detailing_for_abaqus
# Do not run adjust detailing in testing mode as it is time-consuming and is not required to test this function
if not project.rhdhvABQ.test_mode:
_adjust_detailing_for_abaqus(project)
force_output_request = {
'directions': [1, 2, 3],
'items': []}
for fstrip in project.collections.fstrips:
element_set_names = []
coordinates = []
if isinstance(fstrip.contour, Polyline):
connecting_shapes = fstrip.get_connecting_shapes()
coordinates = []
for shape in connecting_shapes:
if isinstance(shape, LineMass):
continue
if shape in project.collections.fstrips:
continue
lines = shape.get_connecting_lines(shape=fstrip, include_openings=False)
if lines:
instance, _set = project.rhdhvABQ.get_part_set_name(shape)
element_set_names.append(f"{instance}.{_set}")
coordinates.append([])
# Sort lines for stable tests
lines.sort(key=lambda x: x.name)
for line in lines:
coordinates[-1].append([line.node_start.coordinates, line.node_end.coordinates])
else:
if isinstance(shape, Lines):
points = shape.get_connecting_nodes(shape=fstrip, include_openings=False)
instance, _set = project.rhdhvABQ.get_part_set_name(shape)
if len(points) == 1:
coordinates.append([[points[0].coordinates]])
element_set_names.append(f"{instance}.{_set}")
else:
continue
location = {
'rhdhv_fem': {
'shape_name': fstrip.name,
'geometry_name:': fstrip.contour.name},
'abaqus': {
'element_set_name': element_set_names,
'geometry_coordinates': coordinates}}
force_output_request['items'].append(location)
displacement_output_request = {
'directions': [1, 2, 3],
'items': []}
for pile in project.collections.piles:
instance, _set = project.rhdhvABQ.get_part_set_name(pile)
nodes = pile.get_nodes()
top_node = None
bottom_node = None
if fem_greater(nodes[0].coordinates[2], nodes[-1].coordinates[2]):
top_node = nodes[0]
bottom_node = nodes[-1]
elif fem_smaller(nodes[0].coordinates[2], nodes[-1].coordinates[2]):
top_node = nodes[-1]
bottom_node = nodes[0]
location_force = {
'rhdhv_fem': {
'shape_name': pile.name,
'geometry_name:': bottom_node.name},
'abaqus': {
'element_set_name': [f'{instance}.{_set}'],
'geometry_coordinates': [[[bottom_node.coordinates]]]}}
force_output_request['items'].append(location_force)
if pile_displacement_requets:
location_displacement = {
'rhdhv_fem': {
'shape_name': pile.name,
'geometry_name:': bottom_node.name},
'abaqus': {
'element_set_name': [f'{instance}.{_set}'],
'geometry_coordinates': [[[top_node.coordinates, bottom_node.coordinates]]]}}
displacement_output_request['items'].append(location_displacement)
with open(json_path, 'r') as f:
j = json.load(f)
with open(json_path, 'w') as f:
j['section_force_requests'] = force_output_request
if pile_displacement_requets:
j['displacement_requests'] = displacement_output_request
json.dump(j, f, indent=4)
### ===================================================================================================================
### 3. End of Script
### ===================================================================================================================