Hi @dokken. Thank you very much for your time.
I found a very similar issue in Error in Result at boundary for a Poisson Subdomains/Interior boundary problem - FEniCS Project
Apparently this might be a bug in assemble_system function in dolfin v2019 - raised here - fenics-project / DOLFIN / issues / #1075 - Assembly of integrals over interal boundaries messes up DirichletBCs — Bitbucket
I had to replace NonLinearVariationalProblem with custom newton solver - Set krylov linear solver paramters in newton solver - #3 by nate
I am still unsure of what would be the safest option when dealing with nonlinear problems with integration over interior boundaries.