|
| 1 | +from __future__ import division |
| 2 | + |
| 3 | +import numpy as np |
| 4 | +import precice |
| 5 | + |
| 6 | + |
| 7 | +def displace_flap(x, y, t, flap_tip_y): |
| 8 | + x_displ = np.zeros_like(x) |
| 9 | + y_displ = np.zeros_like(y) |
| 10 | + # get displacement independent of x, only dependent on y and t |
| 11 | + max_x_displ = 0.5 |
| 12 | + period_fac = 3 * np.pi |
| 13 | + damping_fac = 8 # damps the amplitude of the sine |
| 14 | + # defines how much the sine is shifted in y-direction |
| 15 | + shift = 0.95 |
| 16 | + # wiggles the flap periodically. |
| 17 | + # the arcsin(-shift) in the sine evaluation is necessary to start at a flap displacement of 0 at t=0 |
| 18 | + # (i.e. sin(arcsin(-shift))+shift = 0) |
| 19 | + x_displ = np.minimum(((np.sin(period_fac * t + np.arcsin(-shift)) + shift) / |
| 20 | + damping_fac), max_x_displ) * y / flap_tip_y |
| 21 | + |
| 22 | + dimensions = 2 |
| 23 | + displ = np.zeros((len(x), dimensions)) |
| 24 | + displ[:, 0] = x_displ |
| 25 | + displ[:, 1] = y_displ |
| 26 | + |
| 27 | + return displ |
| 28 | + |
| 29 | + |
| 30 | +configuration_file_name = "../precice-config.xml" |
| 31 | +participant_name = "Solid" |
| 32 | +mesh_name = "Solid-Mesh" |
| 33 | +write_data_name = 'Displacement' |
| 34 | + |
| 35 | +solver_process_index = 0 |
| 36 | +solver_process_size = 1 |
| 37 | + |
| 38 | +# define mesh |
| 39 | +H = 1 |
| 40 | +W = 0.1 |
| 41 | + |
| 42 | +interface = precice.Participant(participant_name, configuration_file_name, solver_process_index, solver_process_size) |
| 43 | +dimensions = interface.get_mesh_dimensions(mesh_name) |
| 44 | +assert (dimensions == 2) |
| 45 | + |
| 46 | +x_left = 0.0 - 0.5 * W # left boundary of the flap |
| 47 | +x_right = 0.5 * W # right boundary of the flap |
| 48 | +y_bottom = 0.0 # bottom of the flap |
| 49 | +y_top = y_bottom + H # top of the flap |
| 50 | + |
| 51 | +n = 24 # Number of vertices per side |
| 52 | +t = 0 |
| 53 | + |
| 54 | +vertices = np.zeros((2 * n, dimensions)) |
| 55 | +# define vertices of flap's left side |
| 56 | +vertices[:n, 1] = np.linspace(y_bottom, y_top, n) |
| 57 | +vertices[:n, 0] = x_left |
| 58 | +# define vertices of flap's right side |
| 59 | +vertices[n:, 1] = np.linspace(y_bottom, y_top, n) |
| 60 | +vertices[n:, 0] = x_right |
| 61 | + |
| 62 | +vertex_ids = interface.set_mesh_vertices(mesh_name, vertices) |
| 63 | + |
| 64 | +if interface.requires_initial_data(): |
| 65 | + # initially, there should be no displacement |
| 66 | + interface.write_data(np.zeros_like(vertices)) |
| 67 | + |
| 68 | +interface.initialize() |
| 69 | +# change if necessary |
| 70 | +solver_dt = np.inf |
| 71 | +# for checkpointing |
| 72 | +t_cp = 0 |
| 73 | + |
| 74 | +while interface.is_coupling_ongoing(): |
| 75 | + if interface.requires_writing_checkpoint(): |
| 76 | + t_cp = t |
| 77 | + |
| 78 | + precice_dt = interface.get_max_time_step_size() |
| 79 | + dt = min([solver_dt, precice_dt]) |
| 80 | + # wiggle the flap |
| 81 | + write_data = displace_flap(vertices[:, 0], vertices[:, 1], t, H) |
| 82 | + |
| 83 | + interface.write_data(mesh_name, write_data_name, vertex_ids, write_data) |
| 84 | + interface.advance(dt) |
| 85 | + |
| 86 | + if interface.requires_reading_checkpoint(): |
| 87 | + t = t_cp |
| 88 | + else: |
| 89 | + # update t |
| 90 | + t += dt |
| 91 | + |
| 92 | +interface.finalize() |
0 commit comments