fabiangabriel/active_flow_control_past_cylinder_using_drl

active_flow_control_past_cylinder_using_drl's Issues

Stabilize control

Hi Fabian,

the next task will be to stabilize the learning process by constraining the angular velocity. There are two means to restrict the angular velocity:

1. limiting the angular velocity applied as the boundary condition
2. limiting the distribution from which omega is sampled during training

The first point is the most essential one since extreme values of omega ultimately make the simulations crash. Restricting the exploration/the distribution should lead to a more efficient sampling, and therefore, to a smoother and faster training.

To implement step 1., you need to make modifications in the agentRotatingWallVelocityFvPatchVectorField.C file and the get_predictions method of the policy network for computing a normalized omega. You find a reference implementation on page 394 of grokking DRL:

mean, log_std = self.forward(state)
pi_s = Normal(mean, log_std.exp())
pre_tanh_action = pi_s.rsample()
tanh_action = torch.tanh(pre_tanh_action)
action = self.rescale_fn(tanh_action)
log_prob = pi_s.log_prob(pre_tanh_action) - torch.log(
            (1 - tanh_action.pow(2)).clamp(0, 1) + epsilon)
log_prob = log_prob.sum(dim=1, keepdim=True)

Note that the reason for having to implement the same code both in Python and C++ is that, at the time of writing, tracing probability distributions is not supported (read more).

Tip: Before running a training, make sure that both implementations consistently deliver outputs in the prescribed range.

As always, if questions or problems come up, feel free to post them here.

Best, Andre

Create new *blockMeshDict*

Hi Fabian,
as a warm-up, I would like to create an alternative mesh for the cylinder flow setup. The current setup uses a combination of the native OpenFOAM meshers blockMesh and snappyHexMesh. To simplify and accelerate the mesh creation, I would like to create the mesh entirely in blockMesh. I suggest the following steps

• get familiar with the current simulation setup
• check out the OpenFOAM user guide of blockMesh
• this tutorial from the tutorial guide might be helpful
• modify this existing blockMeshDict and remove snappyHexMeshDict and extrudeMeshDict

If questions come up, you can post them in this issue.

Best, Andre

Parallelize agentRotatingWallVelocity boundary condition

Hi Fabian,
in the next step, we allow the cylinder patch to be distributed across processors. You find all required operations in this tutorial. Communication in OpenFOAM is implemented in the PStream class. What you have to do on a top-level:

1. gather list of pressure values on the master process
2. run forward pass through the policy network to obtain omega (on the master only)
3. communicate the new omega to the other processors

The pressure values on the patch are accessed with the statement:

const fvPatchField<scalar> &p = patch().lookupPatchField<volScalarField, scalar>("p");

In this post you find a potential solution to gather values on the master.
The part of the code that has to be executed on the master only is:

            torch::Tensor features = torch::zeros({1, p.size()}, torch::kFloat64);
            forAll(p, i)
            {
                features[0][i] = p[i];
            }
            std::vector<torch::jit::IValue> policyFeatures{features};
            torch::Tensor gauss_parameters = policy_.forward(policyFeatures).toTensor();
            torch::Tensor mean_tensor = gauss_parameters[0][0];
            std::cout << "log_std: " << gauss_parameters[0][1].item<double>() << "\n";
            torch::Tensor log_std_tensor = torch::clamp(gauss_parameters[0][1], -20.0, log_std_max_);
            std::cout << "clipped log_std: " << log_std_tensor << "\n";
            if (train_)
            {
                // sample from Gaussian distribution during training
                omega_ = at::normal(mean_tensor, log_std_tensor.exp()).item<double>();
            }
            else
            {
                // use expected (mean) angular velocity
                omega_ = mean_tensor.item<double>();
            }

Instead of p.size() you have to first determine the overall number of faces by summing up the individual number of faces on each processor (the master has a list of lists; the master sums up the lengths of all lists). Once you know the number of faces, you can allocate the features tensor, assign the values from all lists to the tensor, and run the forward pass. After omega_ was assigned on the master, you still have to scatter it to the other processors. Finally, also the code

scalar mean = mean_tensor.item<double>();
scalar log_std = log_std_tensor.item<double>();
scalar var = (log_std_tensor + log_std_tensor).exp().item<double>();
scalar entropy = 0.5 + 0.5*log(2.0*M_PI) + log_std;
scalar log_p = -((omega_ - mean) * (omega_ - mean)) / (2.0*var) - log_std - log(sqrt(2.0*M_PI));
saveTrajectory(log_p, entropy, mean, log_std);
// reset cumulative values
theta_cumulative_ = 0.0;
dt_theta_cumulative_ = 0.0;
Info << "New omega: " << omega_ << "; old value: " << omega_old_ << "\n";

has to be executed by the master only. The line update_omega_ = false; has to be executed by all processes.
If I missed something or problems come up, let me know.

Edit: I just realized that also the saveTrajectory function needs some modification since it writes out all pressure values. I suggest passing the features tensor to the function since it already contains the gathered pressure values.

Best, Andre

Training does shutdown right after start

Hey Fabian, I got a problem when I execute the python_job.sh in order to start a training:
After I execute the script I can see the 12 trajectories are being generated in the JOB screen but they are closed right after start and the SLURM output files just say this here:

/var/tmp/slurmd_spool/job1711768/slurm_script: Zeile 13: ./Allrun.singularity: Keine Berechtigung

Do you know any fix for this?
Maybe I set up the repository incorrectly..

greetings
Erik

Mesh dependency study

Hi Fabian,

I suggest a mesh dependency study for different Reynolds numbers as a next step to advance the project. As you start using denser meshes, the computational demand increases. You can accelerate the simulations by increasing the number of processes. I would impose an upper limit of 10 processes per simulation. To change the number of processes, you need to modify three files:

1. Allrun and Allrun.singularity
2. system/decomposeParDict

As an upper limited for the Reynolds number, I suggest Re=1000. To change the Reynolds number, modify the inflow boundary condition in 0.org/U. I would conduct the mesh debendency study roughly as follows:

• create three or four different mesh densities; each new density level should have twice as many cells per direction as the previous one; use 10 processes by default
• set up cases for several different Reynolds numbers between 100 and 1000 and execute them; the Courant number should stay roughly constant, so the timestep in system/controlDict needs adjustment depending on the inflow velocity and the mesh density; you can see the Courant number in log.pimpleFoam
• analyze how lift and drag change with varying mesh density; this notebook might be a good reference for the visualization
• observe how the mesh dependency increases with increasing Reynolds number.

Best, Andre

Test and copy Darshan's instructions

Hi Fabian,

Darshan has now taken the time to clean up his repository nicely. I would like you to have a look at his instructions in the README file and to see if those instructions also work for you. If you encounter any trouble or if there are missing links, which is likely going to happen, please open an issue in Darshan's repository such that the two of you can find a solution together (or the three of us if needed). Once you can reproduce Darshan's DRL workflow, I would like you to copy the DRL-specific parts of the repo, code and instructions, to your repository. Of course, you should add a reference to the origin in the header of each source file. So, in summary, there are two top-level tasks:

• Test Darshan's setup and instructions
• copy the DRL part of the repo and add a reference to the original source

This issue is going to take some time to complete, but it is extremely important. Take your time to check things carefully :-)

Best, Andre