{"action":"created","issue":{"url":"https://api.github.com/repos/thetisproject/thetis/issues/26","repository_url":"https://api.github.com/repos/thetisproject/thetis","labels_url":"https://api.github.com/repos/thetisproject/thetis/issues/26/labels{/name}","comments_url":"https://api.github.com/repos/thetisproject/thetis/issues/26/comments","events_url":"https://api.github.com/repos/thetisproject/thetis/issues/26/events","html_url":"https://github.com/thetisproject/thetis/issues/26","id":143097222,"number":26,"title":"Redesign equations","user":{"login":"tkarna","id":7586775,"avatar_url":"https://avatars.githubusercontent.com/u/7586775?v=3","gravatar_id":"","url":"https://api.github.com/users/tkarna","html_url":"https://github.com/tkarna","followers_url":"https://api.github.com/users/tkarna/followers","following_url":"https://api.github.com/users/tkarna/following{/other_user}","gists_url":"https://api.github.com/users/tkarna/gists{/gist_id}","starred_url":"https://api.github.com/users/tkarna/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/tkarna/subscriptions","organizations_url":"https://api.github.com/users/tkarna/orgs","repos_url":"https://api.github.com/users/tkarna/repos","events_url":"https://api.github.com/users/tkarna/events{/privacy}","received_events_url":"https://api.github.com/users/tkarna/received_events","type":"User","site_admin":false},"labels":[{"url":"https://api.github.com/repos/thetisproject/thetis/labels/enhancement","name":"enhancement","color":"84b6eb"}],"state":"open","locked":false,"assignee":null,"milestone":null,"comments":3,"created_at":"2016-03-23T22:32:11Z","updated_at":"2016-04-25T07:05:54Z","closed_at":null,"body":"# Redesign equations\r\n\r\nThe `Equation` class as it is implemented now is too rigid for developing custom time integrators. We need a more flexible API where the user can choose which terms are treated implicitly/explicitly/differently.\r\n\r\n## Equations and terms\r\n\r\nEquation class defines all the terms in the equation. It takes a bunch of options and the function space where the solution belongs:\r\n\r\n```python\r\nfunction_space = solver.function_spaces.P1DG\r\noptions = {'use_grad_div_viscosity': True}\r\neq = MyEquation(function_space, **options)\r\n\r\nclass MyEquation(Equation):\r\n    def __init__(self, function_space, ...)\r\n        self.function_space = function_space\r\n        self.test = TestFunction(self.function_space)\r\n        ...\r\n```\r\n\r\nOptions are typically booleans like above.\r\n\r\nThe terms are implemented as methods that return the corresponding form, including interface terms and boundary conditions\r\n\r\n```python\r\nclass MyEquation(Equation):\r\n    ...\r\n    def viscosity(self, fields, bnd_conditions):\r\n        \"\"\"viscosity operator with interior penalty terms.\"\"\"\r\n        grad_div = self.options.use_grad_div_viscosity\r\n        f = 0\r\n        f += ... self.test * dx\r\n        return f\r\n```\r\n\r\nSign convention: all terms are assumed to be on the left hand side of the equation A + term = 0.\r\n\r\nThe terms depend on a `fields` dictionary\r\n\r\n```python\r\nfields={'uv_3d': solver.fields.uv_3d,\r\n        'elev_3d': solver.fields.elev_3d,\r\n        'viscosity': None}\r\n```\r\n\r\nFields should be `Coefficient`s, e.g. `Function` or `Constant` instances, or UFL expressions(?). Some fields are optional, e.g. if `viscosity` is omitted (or it is `None`) the corresponding terms will be ignored. (Alternatively we could have `option.use_viscosity`.) Some fields are required and should be asserted.\r\n\r\nThe equation object knows all the terms that it implements\r\n\r\n```python\r\nclass MyEquation(Equation):\r\n    def __init__(...):\r\n        ...\r\n        # add some python magic to do this\r\n        self.terms = {}\r\n        self.terms['viscosity'] = self.viscosity\r\n        self.terms['elevation_gradient'] = self.elevation_gradient\r\n```\r\n\r\nTo differentiate how the terms should be treated, the user can add labels:\r\n\r\n```python\r\neq.set_label('implicit', ('viscosity', 'elevation_gradient'))\r\neq.set_label('explicit', 'advection')\r\neq.set_label('source', 'wind_stress')\r\n```\r\n\r\nThe equation class provides a method for evaluating the forms\r\n\r\n```\r\nf_implicit = eq.get_form('implicit', fields)\r\n# returns eq.viscosity(fields) + eq.elevation_gradient(fields)\r\n\r\nf_explicit = eq.get_form('explicit', fields)\r\n# returns eq.advection(fields)\r\n```\r\n\r\n## Time integrator\r\n\r\nTime integrator then takes the equation, fields and labels as an argument \r\n\r\n```python\r\nlabels = {'implicit': ('viscosity', 'elevation_gradient'),\r\n          'explicit': 'advection'}\r\nti = MyTimeIntegrator(eq, fields, labels, bnd_conditions)\r\n```\r\n\r\nThe labels depend on the time integrator, e.g. backward Euler only has implicit terms.\r\nUsing the `Equation.get_form` method the time integrator can then construct the solvers it needs:\r\n\r\n```python\r\nclass MyTimeIntegrator(Equation):\r\n    def __init__(self, equation, fields, labels, bnd_conditions):\r\n        self.equation = equation\r\n        ...\r\n        mass_term = ...\r\n        f_k1 = mass_term + self.dt*self.equation.get_form('implicit', fields_k1)\r\n        self.problem_k1 = NonlinearVariationalProblem(f_k1, self.solution_k1)\r\n        self.solver_k1 = NonlinearVariationalSolver(self.problem_k1)\r\n\r\n    def set_timestep(self, dt):\r\n        self.dt.assign(dt)\r\n\r\n    def advance(self, time, update_forcings=None):\r\n        ...\r\n```\r\n        \r\nNOTE: Currently the equation object knows the fields it uses and also the solution function. These however may depend on the time integrator so it's better to define them here.\r\n\r\n## Boundary conditions\r\n\r\nBoundary conditions are defined in a dictionary that define which fields should be replaced at the boundary (for weak boundary conditions). There can be additional boundary specific options.\r\n\r\n\r\n```python\r\nocean_bnd_id = 1\r\nriver_bnd_id = 2\r\nocean_uv = Constant((1, 0, 0))\r\nriver_uv = Constant((0.05, 0, 0))\r\nocean_funcs = {'uv_3d': ocean_uv}\r\nriver_funcs = {'uv_3d': river_uv, 'clamped': True}\r\nbnd_advection = {ocean_bnd_id: ocean_funcs,\r\n                 river_bnd_id: river_funcs}\r\n    \r\n# apply to term\r\nf_adv = eq.advection(fields, bnd_advection)\r\n\r\n# or time integrator\r\nbcs = {'advection': bnd_advection}\r\nti = MyTimeIntegrator(eq, fields, labels, bnd_conditions)\r\n```\r\n\r\nIt should be possible to support Dirichlet boundary conditions in a similar fashion,\r\n```python\r\nocean_funcs = DirichletBC(...)\r\n```\r\nthis just will be passed to the solver instead of the term.\r\n\r\n## Naming convention\r\n\r\nFor the sake of clarity the field names should follow `field_metadata` (defined in [field_defs.py](../blob/master/thetis/field_defs.py) ) as closely as possible. Likewise, options should be similar to [options.py](../blob/master/thetis/options.py) where applicable.\r\n"},"comment":{"url":"https://api.github.com/repos/thetisproject/thetis/issues/comments/214172346","html_url":"https://github.com/thetisproject/thetis/issues/26#issuecomment-214172346","issue_url":"https://api.github.com/repos/thetisproject/thetis/issues/26","id":214172346,"user":{"login":"tkarna","id":7586775,"avatar_url":"https://avatars.githubusercontent.com/u/7586775?v=3","gravatar_id":"","url":"https://api.github.com/users/tkarna","html_url":"https://github.com/tkarna","followers_url":"https://api.github.com/users/tkarna/followers","following_url":"https://api.github.com/users/tkarna/following{/other_user}","gists_url":"https://api.github.com/users/tkarna/gists{/gist_id}","starred_url":"https://api.github.com/users/tkarna/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/tkarna/subscriptions","organizations_url":"https://api.github.com/users/tkarna/orgs","repos_url":"https://api.github.com/users/tkarna/repos","events_url":"https://api.github.com/users/tkarna/events{/privacy}","received_events_url":"https://api.github.com/users/tkarna/received_events","type":"User","site_admin":false},"created_at":"2016-04-25T07:05:54Z","updated_at":"2016-04-25T07:05:54Z","body":"Implemented in #35 "}}