# Path: dolfin_navier_scipy/residual_checks.py
# def get_steady_state_res(V=None, outflowds=None, gradvsymmtrc=True, nu=None):
# 
#     def steady_state_res(vel, pres, phi=None):
#         if phi is None:
#             phi = dolfin.TestFunction(V)
# 
#         cnvfrm = inner(dolfin.dot(vel, nabla_grad(vel)), phi)*dx
#         diffrm = nu*inner(grad(vel)+grad(vel).T, grad(phi))*dx
#         if gradvsymmtrc:
#             nvec = dolfin.FacetNormal(V.mesh())
#             diffrm = diffrm - (nu*inner(grad(vel).T*nvec, phi))*outflowds
# 
#         pfrm = inner(pres, div(phi))*dx
#         res = dolfin.assemble(diffrm+cnvfrm-pfrm)
#         return res
# 
#     return steady_state_res


import numpy as np
import scipy.optimize as sco
import dolfin
import dolfin_navier_scipy.stokes_navier_utils as snu
import dolfin_navier_scipy.dolfin_to_sparrays as dts
import dolfin_navier_scipy.problem_setups as dnsps
import dolfin_navier_scipy.data_output_utils as dou
from dolfin_navier_scipy.residual_checks import get_steady_state_res
                   dbcinds=femp['dbcinds'], dbcvals=femp['dbcvals'],
                   data_prfx=ddir+data_prfx,
                   paraviewoutput=ParaviewOutput,
                   vfileprfx=proutdir+'vel_',
                   pfileprfx=proutdir+'p_')
    L = femp['charlen']  # characteristic length

    phionevec = np.zeros((femp['V'].dim(), 1))
    phionevec[femp['mvwbcinds'], :] = 1.
    phione = dolfin.Function(femp['V'])
    phione.vector().set_local(phionevec)
    pickx = dolfin.as_matrix([[1., 0.], [0., 0.]])
    picky = dolfin.as_matrix([[0., 0.], [0., 1.]])
    pox = pickx*phione
    poy = picky*phione

    phitwovec = np.zeros((femp['V'].dim(), 1))
    phitwovec[femp['mvwbcinds'], 0] = femp['mvwbcvals']
    phitwo = dolfin.Function(femp['V'])
    phitwo.vector().set_local(phitwovec)
    if ParaviewOutput:
        phifile = dolfin.File('results/phione.pvd')
        phifile << phitwo

    # getld = dnsps.LiftDragSurfForce(V=femp['V'], nu=nu,
    #                                 phione=phione, phitwo=phitwo,
    #                                 outflowds=femp['outflowds'],
    #                                 ldds=femp['liftdragds'])

    steady_state_res = \

# Path: dolfin_navier_scipy/residual_checks.py
# def get_steady_state_res(V=None, outflowds=None, gradvsymmtrc=True, nu=None):
# 
#     def steady_state_res(vel, pres, phi=None):
#         if phi is None:
#             phi = dolfin.TestFunction(V)
# 
#         cnvfrm = inner(dolfin.dot(vel, nabla_grad(vel)), phi)*dx
#         diffrm = nu*inner(grad(vel)+grad(vel).T, grad(phi))*dx
#         if gradvsymmtrc:
#             nvec = dolfin.FacetNormal(V.mesh())
#             diffrm = diffrm - (nu*inner(grad(vel).T*nvec, phi))*outflowds
# 
#         pfrm = inner(pres, div(phi))*dx
#         res = dolfin.assemble(diffrm+cnvfrm-pfrm)
#         return res
# 
#     return steady_state_res


# Path: tests/steadystate_rotcyl.py
import numpy as np
import scipy.optimize as sco
import dolfin
import dolfin_navier_scipy.stokes_navier_utils as snu
import dolfin_navier_scipy.dolfin_to_sparrays as dts
import dolfin_navier_scipy.problem_setups as dnsps
import dolfin_navier_scipy.data_output_utils as dou
from dolfin_navier_scipy.residual_checks import get_steady_state_res

                   dbcinds=femp['dbcinds'], dbcvals=femp['dbcvals'],
                   data_prfx=ddir+data_prfx,
                   paraviewoutput=ParaviewOutput,
                   vfileprfx=proutdir+'vel_',
                   pfileprfx=proutdir+'p_')
    L = femp['charlen']  # characteristic length

    phionevec = np.zeros((femp['V'].dim(), 1))
    phionevec[femp['mvwbcinds'], :] = 1.
    phione = dolfin.Function(femp['V'])
    phione.vector().set_local(phionevec)
    pickx = dolfin.as_matrix([[1., 0.], [0., 0.]])
    picky = dolfin.as_matrix([[0., 0.], [0., 1.]])
    pox = pickx*phione
    poy = picky*phione

    phitwovec = np.zeros((femp['V'].dim(), 1))
    phitwovec[femp['mvwbcinds'], 0] = femp['mvwbcvals']
    phitwo = dolfin.Function(femp['V'])
    phitwo.vector().set_local(phitwovec)
    if ParaviewOutput:
        phifile = dolfin.File('results/phione.pvd')
        phifile << phitwo

    # getld = dnsps.LiftDragSurfForce(V=femp['V'], nu=nu,
    #                                 phione=phione, phitwo=phitwo,
    #                                 outflowds=femp['outflowds'],
    #                                 ldds=femp['liftdragds'])

    steady_state_res = \

        get_steady_state_res(V=femp['V'], gradvsymmtrc=True,