{"action":"created","comment":{"url":"https://api.github.com/repos/geodynamics/aspect/pulls/comments/30956849","id":30956849,"diff_hunk":"@@ -0,0 +1,444 @@\n+/*\n+  Copyright (C) 2011 - 2015 by the authors of the ASPECT code.\n+\n+  This file is part of ASPECT.\n+\n+  ASPECT is free software; you can redistribute it and/or modify\n+  it under the terms of the GNU General Public License as published by\n+  the Free Software Foundation; either version 2, or (at your option)\n+  any later version.\n+\n+  ASPECT is distributed in the hope that it will be useful,\n+  but WITHOUT ANY WARRANTY; without even the implied warranty of\n+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the\n+  GNU General Public License for more details.\n+\n+  You should have received a copy of the GNU General Public License\n+  along with ASPECT; see the file doc/COPYING.  If not see\n+  <http://www.gnu.org/licenses/>.\n+*/\n+\n+#include <aspect/material_model/depth_dependent.h>\n+#include <deal.II/base/parameter_handler.h>\n+#include <utility>\n+#include <limits>\n+using namespace dealii;\n+\n+namespace aspect\n+{\n+  namespace MaterialModel\n+  {\n+\n+    namespace internal\n+    {\n+      class RadialViscosityLookup\n+      {\n+        public:\n+          RadialViscosityLookup(const std::string &filename)\n+          {\n+            /* This method is used for the Table method of depth dependent viscosity */\n+            std::string temp;\n+            std::ifstream in(filename.c_str(), std::ios::in);\n+            AssertThrow (in,\n+                         ExcMessage (std::string(\"Couldn't open file <\") + filename));\n+\n+            min_depth=DBL_MAX;\n+            max_depth=DBL_MIN;\n+\n+            // The input viscosity file has two columns that are the viscosity and Depth\n+            std::string header;\n+            getline(in, header);/* Discard header line */\n+            while (!in.eof())\n+              {\n+                double visc, depth;\n+                in >> visc;;\n+                if (in.eof())\n+                  break;\n+                in >> depth;\n+                getline(in, temp);\n+\n+                min_depth = std::min(depth, min_depth);\n+                max_depth = std::max(depth, max_depth);\n+\n+                values.push_back(visc);\n+                depthvalues.push_back(depth);\n+              }\n+          }\n+\n+          double radial_viscosity(double depth)\n+          {\n+            // Interpolate the viscosity to the nearest neighbour\n+            std::vector<double> depth_diff (values.size(), 0);\n+            Assert(depth>=min_depth, ExcMessage(\"not in range\"));\n+            Assert(depth<=max_depth, ExcMessage(\"not in range\"));\n+\n+            for (int i = 0; i < values.size(); i++)\n+              depth_diff[i] = std::abs(depthvalues[i] - depth);\n+\n+            double depth_val = 1e7;\n+            for (int i = 0; i < values.size(); i++)\n+              depth_val = std::min(depth_diff[i],depth_val);\n+\n+            unsigned int idx = values.size();\n+            for (int i = 0; i < values.size(); i++)\n+              if (depth_val == std::abs(depthvalues[i] - depth))\n+                idx = i;\n+\n+            Assert(idx<values.size(), ExcMessage(\"not in range\"));\n+            return values[idx];\n+          }\n+\n+        private:\n+          std::vector<double> depthvalues;\n+          std::vector<double> values;\n+          double min_depth;\n+          double max_depth;\n+      };\n+    }\n+\n+    template <int dim>\n+    double\n+    DepthDependent<dim>::viscosity_from_list(const double &depth) const\n+    {\n+      /* find the layer containing the specified depth and return the corresponding viscosity */\n+      int i=0;\n+      const int nlayers = depth_values.size()-1;\n+      while ( depth > depth_values[i] && i<nlayers ) i++; /* incerment i until depth is above base of layer i */\n+      return viscosity_values[i];\n+    }\n+\n+    template <int dim>\n+    void\n+    DepthDependent<dim>::initialize()\n+    {\n+      /* If using the File method for depth-dependence, initialize the lookup table */\n+      if ( viscosity_method == File )\n+        {\n+          radial_viscosity_lookup.reset(new internal::RadialViscosityLookup(datadirectory+radial_viscosity_file_name));\n+        }\n+    }\n+\n+    template <int dim>\n+    double\n+    DepthDependent<dim>::calculate_depth_dependent_prefactor(const double &depth) const\n+    {\n+      /* The depth dependent prefactor is the multiplicative factor by which the\n+       * viscosity computed by the base model's evaluate mathod will be scaled */\n+      const double reference_viscosity = base_model->reference_viscosity();\n+      if ( viscosity_method == File )\n+        {\n+          return radial_viscosity_lookup->radial_viscosity(depth)/reference_viscosity;\n+        }\n+      else if ( viscosity_method == Function )\n+        {\n+          Point<1> dpoint(depth);\n+          const double viscosity_depth_prefactor = viscosity_function.value(dpoint);\n+          Assert (viscosity_depth_prefactor > 0.0, ExcMessage(\"Viscosity depth function should be larger than zero\"));\n+          return viscosity_depth_prefactor/reference_viscosity;\n+        }\n+      else if (viscosity_method == List)\n+        {\n+          return viscosity_from_list(depth)/reference_viscosity;\n+        }\n+      else if (viscosity_method == None)\n+        {\n+          return 1.0;\n+        }\n+      else\n+        {\n+          Assert( false, ExcMessage(\"Invalid method for viscosity depth dependence\") );\n+          return 0.0;\n+        }\n+    }\n+\n+\n+    template <int dim>\n+    void\n+    DepthDependent<dim>::evaluate(const typename Interface<dim>::MaterialModelInputs &in,\n+                                  typename Interface<dim>::MaterialModelOutputs &out) const\n+    {\n+      base_model -> evaluate(in,out);\n+      /* Scale the base model viscosity value by the depth dependenct prefactor */\n+      for (unsigned int i=0; i < out.viscosities.size(); ++i)\n+        {\n+          const double depth = this->get_geometry_model().depth(in.position[i]);\n+          out.viscosities[i] *= calculate_depth_dependent_prefactor( depth );\n+        }\n+    }\n+\n+\n+    template <int dim>\n+    void\n+    DepthDependent<dim>::declare_parameters (ParameterHandler &prm)\n+    {\n+      prm.enter_subsection(\"Material model\");\n+      {\n+        prm.enter_subsection(\"Depth dependent model\");\n+        {\n+          prm.declare_entry(\"Base model\",\"simple\",\n+                            Patterns::Selection(\"simple|simpler|simple compressible|Steinberger|table|composition reaction|Morency and Doin|latent heat|latent heat melt|multicomponent\"),","path":"source/material_model/depth_dependent.cc","position":179,"original_position":179,"commit_id":"348904e329cb525e96614da5d9ff4c145c4ee332","original_commit_id":"348904e329cb525e96614da5d9ff4c145c4ee332","user":{"login":"maxrudolph","id":8154934,"avatar_url":"https://avatars.githubusercontent.com/u/8154934?v=3","gravatar_id":"","url":"https://api.github.com/users/maxrudolph","html_url":"https://github.com/maxrudolph","followers_url":"https://api.github.com/users/maxrudolph/followers","following_url":"https://api.github.com/users/maxrudolph/following{/other_user}","gists_url":"https://api.github.com/users/maxrudolph/gists{/gist_id}","starred_url":"https://api.github.com/users/maxrudolph/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/maxrudolph/subscriptions","organizations_url":"https://api.github.com/users/maxrudolph/orgs","repos_url":"https://api.github.com/users/maxrudolph/repos","events_url":"https://api.github.com/users/maxrudolph/events{/privacy}","received_events_url":"https://api.github.com/users/maxrudolph/received_events","type":"User","site_admin":false},"body":"@bangerth Could you please tell me the name of this function?","created_at":"2015-05-24T23:12:36Z","updated_at":"2015-05-24T23:12:36Z","html_url":"https://github.com/geodynamics/aspect/pull/428#discussion_r30956849","pull_request_url":"https://api.github.com/repos/geodynamics/aspect/pulls/428","_links":{"self":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/comments/30956849"},"html":{"href":"https://github.com/geodynamics/aspect/pull/428#discussion_r30956849"},"pull_request":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/428"}}},"pull_request":{"url":"https://api.github.com/repos/geodynamics/aspect/pulls/428","id":36109336,"html_url":"https://github.com/geodynamics/aspect/pull/428","diff_url":"https://github.com/geodynamics/aspect/pull/428.diff","patch_url":"https://github.com/geodynamics/aspect/pull/428.patch","issue_url":"https://api.github.com/repos/geodynamics/aspect/issues/428","number":428,"state":"open","locked":false,"title":"Add depth-dependent viscosity","user":{"login":"maxrudolph","id":8154934,"avatar_url":"https://avatars.githubusercontent.com/u/8154934?v=3","gravatar_id":"","url":"https://api.github.com/users/maxrudolph","html_url":"https://github.com/maxrudolph","followers_url":"https://api.github.com/users/maxrudolph/followers","following_url":"https://api.github.com/users/maxrudolph/following{/other_user}","gists_url":"https://api.github.com/users/maxrudolph/gists{/gist_id}","starred_url":"https://api.github.com/users/maxrudolph/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/maxrudolph/subscriptions","organizations_url":"https://api.github.com/users/maxrudolph/orgs","repos_url":"https://api.github.com/users/maxrudolph/repos","events_url":"https://api.github.com/users/maxrudolph/events{/privacy}","received_events_url":"https://api.github.com/users/maxrudolph/received_events","type":"User","site_admin":false},"body":"from Jacqueline Austermann (@jaustermann) and Max Rudolph at ASPECT Hackathon 2015\r\n\r\nWe add a new material model that implements depth-dependent viscosity. Our material model includes documentation, an example data file containing a depth-dependent viscosity structure, and tests of three different methods by which the depth dependence can be specified in 3D spherical shell geometry as well as 2D cartesian geometry.","created_at":"2015-05-24T22:16:59Z","updated_at":"2015-05-24T23:12:36Z","closed_at":null,"merged_at":null,"merge_commit_sha":"8b54ca5b79db38049fb5267823ecb9c0b9342f0e","assignee":null,"milestone":null,"commits_url":"https://api.github.com/repos/geodynamics/aspect/pulls/428/commits","review_comments_url":"https://api.github.com/repos/geodynamics/aspect/pulls/428/comments","review_comment_url":"https://api.github.com/repos/geodynamics/aspect/pulls/comments{/number}","comments_url":"https://api.github.com/repos/geodynamics/aspect/issues/428/comments","statuses_url":"https://api.github.com/repos/geodynamics/aspect/statuses/348904e329cb525e96614da5d9ff4c145c4ee332","head":{"label":"maxrudolph:depth_dependent_3","ref":"depth_dependent_3","sha":"348904e329cb525e96614da5d9ff4c145c4ee332","user":{"login":"maxrudolph","id":8154934,"avatar_url":"https://avatars.githubusercontent.com/u/8154934?v=3","gravatar_id":"","url":"https://api.github.com/users/maxrudolph","html_url":"https://github.com/maxrudolph","followers_url":"https://api.github.com/users/maxrudolph/followers","following_url":"https://api.github.com/users/maxrudolph/following{/other_user}","gists_url":"https://api.github.com/users/maxrudolph/gists{/gist_id}","starred_url":"https://api.github.com/users/maxrudolph/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/maxrudolph/subscriptions","organizations_url":"https://api.github.com/users/maxrudolph/orgs","repos_url":"https://api.github.com/users/maxrudolph/repos","events_url":"https://api.github.com/users/maxrudolph/events{/privacy}","received_events_url":"https://api.github.com/users/maxrudolph/received_events","type":"User","site_admin":false},"repo":{"id":35979322,"name":"aspect","full_name":"maxrudolph/aspect","owner":{"login":"maxrudolph","id":8154934,"avatar_url":"https://avatars.githubusercontent.com/u/8154934?v=3","gravatar_id":"","url":"https://api.github.com/users/maxrudolph","html_url":"https://github.com/maxrudolph","followers_url":"https://api.github.com/users/maxrudolph/followers","following_url":"https://api.github.com/users/maxrudolph/following{/other_user}","gists_url":"https://api.github.com/users/maxrudolph/gists{/gist_id}","starred_url":"https://api.github.com/users/maxrudolph/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/maxrudolph/subscriptions","organizations_url":"https://api.github.com/users/maxrudolph/orgs","repos_url":"https://api.github.com/users/maxrudolph/repos","events_url":"https://api.github.com/users/maxrudolph/events{/privacy}","received_events_url":"https://api.github.com/users/maxrudolph/received_events","type":"User","site_admin":false},"private":false,"html_url":"https://github.com/maxrudolph/aspect","description":"A parallel, extendible finite element code to simulate thermal convection in both 2D and 3D models.","fork":true,"url":"https://api.github.com/repos/maxrudolph/aspect","forks_url":"https://api.github.com/repos/maxrudolph/aspect/forks","keys_url":"https://api.github.com/repos/maxrudolph/aspect/keys{/key_id}","collaborators_url":"https://api.github.com/repos/maxrudolph/aspect/collaborators{/collaborator}","teams_url":"https://api.github.com/repos/maxrudolph/aspect/teams","hooks_url":"https://api.github.com/repos/maxrudolph/aspect/hooks","issue_events_url":"https://api.github.com/repos/maxrudolph/aspect/issues/events{/number}","events_url":"https://api.github.com/repos/maxrudolph/aspect/events","assignees_url":"https://api.github.com/repos/maxrudolph/aspect/assignees{/user}","branches_url":"https://api.github.com/repos/maxrudolph/aspect/branches{/branch}","tags_url":"https://api.github.com/repos/maxrudolph/aspect/tags","blobs_url":"https://api.github.com/repos/maxrudolph/aspect/git/blobs{/sha}","git_tags_url":"https://api.github.com/repos/maxrudolph/aspect/git/tags{/sha}","git_refs_url":"https://api.github.com/repos/maxrudolph/aspect/git/refs{/sha}","trees_url":"https://api.github.com/repos/maxrudolph/aspect/git/trees{/sha}","statuses_url":"https://api.github.com/repos/maxrudolph/aspect/statuses/{sha}","languages_url":"https://api.github.com/repos/maxrudolph/aspect/languages","stargazers_url":"https://api.github.com/repos/maxrudolph/aspect/stargazers","contributors_url":"https://api.github.com/repos/maxrudolph/aspect/contributors","subscribers_url":"https://api.github.com/repos/maxrudolph/aspect/subscribers","subscription_url":"https://api.github.com/repos/maxrudolph/aspect/subscription","commits_url":"https://api.github.com/repos/maxrudolph/aspect/commits{/sha}","git_commits_url":"https://api.github.com/repos/maxrudolph/aspect/git/commits{/sha}","comments_url":"https://api.github.com/repos/maxrudolph/aspect/comments{/number}","issue_comment_url":"https://api.github.com/repos/maxrudolph/aspect/issues/comments{/number}","contents_url":"https://api.github.com/repos/maxrudolph/aspect/contents/{+path}","compare_url":"https://api.github.com/repos/maxrudolph/aspect/compare/{base}...{head}","merges_url":"https://api.github.com/repos/maxrudolph/aspect/merges","archive_url":"https://api.github.com/repos/maxrudolph/aspect/{archive_format}{/ref}","downloads_url":"https://api.github.com/repos/maxrudolph/aspect/downloads","issues_url":"https://api.github.com/repos/maxrudolph/aspect/issues{/number}","pulls_url":"https://api.github.com/repos/maxrudolph/aspect/pulls{/number}","milestones_url":"https://api.github.com/repos/maxrudolph/aspect/milestones{/number}","notifications_url":"https://api.github.com/repos/maxrudolph/aspect/notifications{?since,all,participating}","labels_url":"https://api.github.com/repos/maxrudolph/aspect/labels{/name}","releases_url":"https://api.github.com/repos/maxrudolph/aspect/releases{/id}","created_at":"2015-05-20T23:18:40Z","updated_at":"2015-05-20T23:18:46Z","pushed_at":"2015-05-24T22:14:19Z","git_url":"git://github.com/maxrudolph/aspect.git","ssh_url":"git@github.com:maxrudolph/aspect.git","clone_url":"https://github.com/maxrudolph/aspect.git","svn_url":"https://github.com/maxrudolph/aspect","homepage":"http://aspect.dealii.org/","size":108770,"stargazers_count":0,"watchers_count":0,"language":"Gnuplot","has_issues":false,"has_downloads":true,"has_wiki":true,"has_pages":false,"forks_count":0,"mirror_url":null,"open_issues_count":0,"forks":0,"open_issues":0,"watchers":0,"default_branch":"master"}},"base":{"label":"geodynamics:master","ref":"master","sha":"8dadd8e8b1d2f5fc7253fcf869c4e5b26778ebec","user":{"login":"geodynamics","id":4052461,"avatar_url":"https://avatars.githubusercontent.com/u/4052461?v=3","gravatar_id":"","url":"https://api.github.com/users/geodynamics","html_url":"https://github.com/geodynamics","followers_url":"https://api.github.com/users/geodynamics/followers","following_url":"https://api.github.com/users/geodynamics/following{/other_user}","gists_url":"https://api.github.com/users/geodynamics/gists{/gist_id}","starred_url":"https://api.github.com/users/geodynamics/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/geodynamics/subscriptions","organizations_url":"https://api.github.com/users/geodynamics/orgs","repos_url":"https://api.github.com/users/geodynamics/repos","events_url":"https://api.github.com/users/geodynamics/events{/privacy}","received_events_url":"https://api.github.com/users/geodynamics/received_events","type":"Organization","site_admin":false},"repo":{"id":19502452,"name":"aspect","full_name":"geodynamics/aspect","owner":{"login":"geodynamics","id":4052461,"avatar_url":"https://avatars.githubusercontent.com/u/4052461?v=3","gravatar_id":"","url":"https://api.github.com/users/geodynamics","html_url":"https://github.com/geodynamics","followers_url":"https://api.github.com/users/geodynamics/followers","following_url":"https://api.github.com/users/geodynamics/following{/other_user}","gists_url":"https://api.github.com/users/geodynamics/gists{/gist_id}","starred_url":"https://api.github.com/users/geodynamics/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/geodynamics/subscriptions","organizations_url":"https://api.github.com/users/geodynamics/orgs","repos_url":"https://api.github.com/users/geodynamics/repos","events_url":"https://api.github.com/users/geodynamics/events{/privacy}","received_events_url":"https://api.github.com/users/geodynamics/received_events","type":"Organization","site_admin":false},"private":false,"html_url":"https://github.com/geodynamics/aspect","description":"A parallel, extendible finite element code to simulate thermal convection in both 2D and 3D models.","fork":false,"url":"https://api.github.com/repos/geodynamics/aspect","forks_url":"https://api.github.com/repos/geodynamics/aspect/forks","keys_url":"https://api.github.com/repos/geodynamics/aspect/keys{/key_id}","collaborators_url":"https://api.github.com/repos/geodynamics/aspect/collaborators{/collaborator}","teams_url":"https://api.github.com/repos/geodynamics/aspect/teams","hooks_url":"https://api.github.com/repos/geodynamics/aspect/hooks","issue_events_url":"https://api.github.com/repos/geodynamics/aspect/issues/events{/number}","events_url":"https://api.github.com/repos/geodynamics/aspect/events","assignees_url":"https://api.github.com/repos/geodynamics/aspect/assignees{/user}","branches_url":"https://api.github.com/repos/geodynamics/aspect/branches{/branch}","tags_url":"https://api.github.com/repos/geodynamics/aspect/tags","blobs_url":"https://api.github.com/repos/geodynamics/aspect/git/blobs{/sha}","git_tags_url":"https://api.github.com/repos/geodynamics/aspect/git/tags{/sha}","git_refs_url":"https://api.github.com/repos/geodynamics/aspect/git/refs{/sha}","trees_url":"https://api.github.com/repos/geodynamics/aspect/git/trees{/sha}","statuses_url":"https://api.github.com/repos/geodynamics/aspect/statuses/{sha}","languages_url":"https://api.github.com/repos/geodynamics/aspect/languages","stargazers_url":"https://api.github.com/repos/geodynamics/aspect/stargazers","contributors_url":"https://api.github.com/repos/geodynamics/aspect/contributors","subscribers_url":"https://api.github.com/repos/geodynamics/aspect/subscribers","subscription_url":"https://api.github.com/repos/geodynamics/aspect/subscription","commits_url":"https://api.github.com/repos/geodynamics/aspect/commits{/sha}","git_commits_url":"https://api.github.com/repos/geodynamics/aspect/git/commits{/sha}","comments_url":"https://api.github.com/repos/geodynamics/aspect/comments{/number}","issue_comment_url":"https://api.github.com/repos/geodynamics/aspect/issues/comments{/number}","contents_url":"https://api.github.com/repos/geodynamics/aspect/contents/{+path}","compare_url":"https://api.github.com/repos/geodynamics/aspect/compare/{base}...{head}","merges_url":"https://api.github.com/repos/geodynamics/aspect/merges","archive_url":"https://api.github.com/repos/geodynamics/aspect/{archive_format}{/ref}","downloads_url":"https://api.github.com/repos/geodynamics/aspect/downloads","issues_url":"https://api.github.com/repos/geodynamics/aspect/issues{/number}","pulls_url":"https://api.github.com/repos/geodynamics/aspect/pulls{/number}","milestones_url":"https://api.github.com/repos/geodynamics/aspect/milestones{/number}","notifications_url":"https://api.github.com/repos/geodynamics/aspect/notifications{?since,all,participating}","labels_url":"https://api.github.com/repos/geodynamics/aspect/labels{/name}","releases_url":"https://api.github.com/repos/geodynamics/aspect/releases{/id}","created_at":"2014-05-06T16:47:18Z","updated_at":"2015-05-08T16:13:02Z","pushed_at":"2015-05-24T23:06:19Z","git_url":"git://github.com/geodynamics/aspect.git","ssh_url":"git@github.com:geodynamics/aspect.git","clone_url":"https://github.com/geodynamics/aspect.git","svn_url":"https://github.com/geodynamics/aspect","homepage":"http://aspect.dealii.org/","size":122783,"stargazers_count":10,"watchers_count":10,"language":"Gnuplot","has_issues":true,"has_downloads":true,"has_wiki":true,"has_pages":false,"forks_count":33,"mirror_url":null,"open_issues_count":56,"forks":33,"open_issues":56,"watchers":10,"default_branch":"master"}},"_links":{"self":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/428"},"html":{"href":"https://github.com/geodynamics/aspect/pull/428"},"issue":{"href":"https://api.github.com/repos/geodynamics/aspect/issues/428"},"comments":{"href":"https://api.github.com/repos/geodynamics/aspect/issues/428/comments"},"review_comments":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/428/comments"},"review_comment":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/comments{/number}"},"commits":{"href":"https://api.github.com/repos/geodynamics/aspect/pulls/428/commits"},"statuses":{"href":"https://api.github.com/repos/geodynamics/aspect/statuses/348904e329cb525e96614da5d9ff4c145c4ee332"}}}}