We present a numerical method to model the dynamics of inextensible
biomembranes in a quasi-Newtonian incompressible flow, which better describes
hemorheology in the small vasculature. We consider a level set model for the
fluid-membrane coupling, while the local inextensibility condition is relaxed
by introducing a penalty term. The penalty method is straightforward to
implement from any Navier-Stokes/level set solver and allows substantial
computational savings over a mixed formulation. A standard Galerkin finite
element framework is used with an arbitrarily high order polynomial
approximation for better accuracy in computing the bending force. The PDE
system is solved using a partitioned strongly coupled scheme based on
Crank-Nicolson time integration. Numerical experiments are provided to validate
and assess the main features of the method