<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Part A
Step1:  Create a 2D Cubic network with standard PSD and define the phase as Water and use Standard physics which implements the washburn capillary pressure relation for throat entry pressure.
Step2:  We can check the model by looking at the model dict on the phys object
Step3:  Now set up and run the algorithm choosing the left and right sides of the network for inlets and outlets respectively. Because we did not set up the network with boundary pores with zero volume a little warning is given because the starting saturation for the algorithm is not zero. However, this is fine and because the network is quite large the starting saturation is actually quite close to zero.
Step4:  The algorithm completes very quickly and the invading phase saturation can be plotted versus the applied boundary pressure.
Step5:  As the network is 2D and cubic we can easily plot the invading phase configuration at the different invasion steps
<ASSISTANT_TASK:>
Python Code:
import openpnm as op
import matplotlib.pyplot as plt
import numpy as np
np.random.seed(10)
from ipywidgets import interact, IntSlider
%matplotlib inline
ws = op.Workspace()
ws.settings["loglevel"] = 40
Explanation: Part A: Ordinary Percolation
OpenPNM contains several percolation algorithms which are central to the multiphase models employed by pore networks. The essential idea is to identify pathways for fluid flow through the network using the entry capillary pressure as a threshold for passage between connected pores. The capillary pressure can either be associated to the pores themselves known as site percolation or the connecting throats known as bond percolation or a mixture of both. OpenPNM provides several models for calculating the entry pressure for a given pore or throat and it generally depends on the size of the pore or throat and the wettability to a particular phase characterised by the contact angle. If a pathway through the network connects pores into clusters that contain both an inlet and an outlet then it is deemed to be percolating.
In this example we will demonstrate Ordinary Percolation which is the fastest and simplest algorithm to run. The number of steps involved in the algorithm is equal to the number of points that are specified in the run method. This can either be an integer, in which case the minimum and maximum capillary entry pressures in the network are used as limits and the integer value is used to create that number of intervals between the limits, or an array of specified pressured can be supplied.
The algorithm progresses incrementally from low pressure to high. At each step, clusters of connected pores are found with entry pressures below the current threshold and those that are not already invaded and connected to an inlet are set to be invaded at this pressure. Therefore the process is quasistatic and represents the steady state saturation that would be achieved if the inlet pressure were to be held at that threshold.
First do our imports
End of explanation
N = 100
net = op.network.Cubic(shape=[N, N, 1], spacing=2.5e-5)
geom = op.geometry.StickAndBall(network=net, pores=net.Ps, throats=net.Ts)
water = op.phases.Water(network=net)
phys = op.physics.Standard(network=net, phase=water, geometry=geom)
Explanation: Create a 2D Cubic network with standard PSD and define the phase as Water and use Standard physics which implements the washburn capillary pressure relation for throat entry pressure.
End of explanation
phys.models['throat.entry_pressure']
Explanation: We can check the model by looking at the model dict on the phys object
End of explanation
#NBVAL_IGNORE_OUTPUT
alg = op.algorithms.OrdinaryPercolation(network=net)
alg.setup(phase=water, pore_volume='pore.volume', throat_volume='throat.volume')
alg.set_inlets(pores=net.pores('left'))
alg.set_outlets(pores=net.pores('right'))
alg.run(points=1000)
alg.plot_intrusion_curve()
plt.show()
Explanation: Now set up and run the algorithm choosing the left and right sides of the network for inlets and outlets respectively. Because we did not set up the network with boundary pores with zero volume a little warning is given because the starting saturation for the algorithm is not zero. However, this is fine and because the network is quite large the starting saturation is actually quite close to zero.
End of explanation
data = alg.get_intrusion_data()
mask = np.logical_and(np.asarray(data.Snwp) > 0.0 , np.asarray(data.Snwp) < 1.0)
mask = np.argwhere(mask).flatten()
pressures = np.asarray(data.Pcap)[mask]
Explanation: The algorithm completes very quickly and the invading phase saturation can be plotted versus the applied boundary pressure.
End of explanation
def plot_saturation(step):
    arg = mask[step]
    Pc = np.ceil(data.Pcap[arg])
    sat = np.around(data.Snwp[arg], 3)
    is_perc = alg.is_percolating(Pc)
    pmask = alg['pore.invasion_pressure'] <= Pc
    im = pmask.reshape([N, N])
    fig, ax = plt.subplots(figsize=[5, 5])
    ax.imshow(im, cmap='Blues');
    title = ('Capillary Pressure: '+str(Pc)+' Saturation: '+str(sat)+
             ' Percolating : '+str(is_perc))
    plt.title(title)
    plt.show()
#NBVAL_IGNORE_OUTPUT
perc_thresh = alg.get_percolation_threshold()
thresh_step = np.argwhere(np.asarray(pressures) == perc_thresh)
interact(plot_saturation, step=IntSlider(min=0, max=len(mask)-1, step=1, value=thresh_step));
Explanation: As the network is 2D and cubic we can easily plot the invading phase configuration at the different invasion steps
End of explanation
<END_TASK>