QuickerSim CFD Toolbox for MATLAB is an incompressible flow solver of Navier-Stokes equations, which works in MATLAB with both a free and full version. It can handle both steady-state and transiet fluid flow simulations. In the free version the fluid dynamics solver code is open and can be viewed by the user to study a bunch of heat transfer and fluid mechanics problems in MATLAB.
Below, we present the script which solves a microfluidic fluid mechanics problem in 3D by means of incompressible Navier-Stokes equations in MATLAB. You can also freely download our toolbox from here and you will find a couple of ‘generic’ solvers that can be applied directly by you to standard CFD problems in MATLAB.
We also show the flow field contour plot, pressure field and concentration field generated in QuickerSim CFD Toolbox for MATLAB for this benchmark microfluidic simulation.
</pre>
<pre><code>clc;
clear;
% Read and display mesh
[p,e,t] = importMeshGmsh('microfluidicH3D_coarse.msh');
displayMesh(p,t);
% Convert mesh to second order
[p,e,t,nVnodes,nPnodes,indices] = convertMeshToSecondOrder(p,e,t);
% Define fluid viscosity
nu = 0.1;
% Init solution and convergence criteria
[u, convergence] = initSolution(p,t,[0 0 0],0);
maxres = 1e-3;
maxiter = 25;
% Define inlet velocity profile
vel = [1 0 0];
% Iterate nonlinear terms
for iter = 1:maxiter
% Assemble matrix and right hand side
[NS, F] = assembleNavierStokesMatrix(p,e,t,nu,u(indices.indu),u(indices.indv),u(indices.indw),'nosupg');
% Apply boundary conditions
[NS, F] = imposeCfdBoundaryCondition(p,e,t,NS,F,[77 78],'inlet', [1 0 0]);
[NS, F] = imposeCfdBoundaryCondition(p,e,t,NS,F,[80 81],'wall');
% Compute and plot residuals
[stop, convergence] = computeResiduals(NS,F,u,size(p),convergence,maxres);
plotResiduals(convergence,2);
% Break if solution converged
if(stop)
break;
end
% Solve equations
u = NS\F;
end
% Plot solution
figure(3);
displaySolution(p,t,u(indices.indu),'x-velocity');
% Solve scalar transport problem
k = 0.001;
% Assemble problem matrix
[D,F] = assembleDiffusionMatrix(p,t,k);
D = D + assembleScalarConvectionMatrix(p,t,k,u(indices.indu),u(indices.indv),u(indices.indw),'supgDoublyAsymptotic');
% Apply boundary conditions
[D,F] = imposeScalarBoundaryCondition(p,e,D,F,77,'value',0);
[D,F] = imposeScalarBoundaryCondition(p,e,D,F,78,'value',1);
% Solve system of equations
phi = D\F;
figure(4);
displaySolution(p,t,phi,'Concentration');
figure(5);
pressure = generatePressureData(u,p,t);
displaySolution(p,t,pressure,'Pressure');
More details can be found in our new Tutorial Center.