```Context
-------

I would like to compute the torque that a  (thin) disc applies onto a ring. I.e. I try to understand what is the impact of this outer ring on
the inner disc in the simulation below.

![Mathematica graphics](http://i.stack.imgur.com/DZtTP.png)

For this I would like to *compute the gravitational potentia*l generated by a (razor) thin disc and a ring.
So the abstraction is the following (seen from 9 different angles)

![Mathematica graphics](http://i.stack.imgur.com/jWVoS.png)

Once I know how to compute the potential it becomes straightforward to compute the torque one feature applies on the other.
I want to use FEM here for flexibility when I will account for a more realistic abstraction of the problem.

Attempt
-------

I have defined a domain

dom = ImplicitRegion[0 <= x <= 1 && -1 <= y <= 1, {x, y}];

and the Laplace operator

op = -Laplacian[u[x, y], {x, theta, y}, "Cylindrical"];

I impose the edge condition that the potential should be 1 on the disc

edge = DirichletCondition[u[x, y] == 1, 0 <= x <= 1/2 && y == 0];

When I solve

uD = NDSolveValue[{op == 0, edge},  u, {x, y} \[Element] dom]

I get

StreamPlot[-{D[#, x], D[#, y]} &@uD[x, y] // Evaluate, {x, y} \[Element] dom,AspectRatio->2]
![Mathematica graphics](http://i.stack.imgur.com/uwjvM.png)

> **Problems**:

> (i) The outer box imposes a (box like) symmetry which is not in the sought solution

> (ii) Strangely enough The code fails if I use the ring-like boundary condition instead:

{DirichletCondition[u[x, y] == 1, 1/2 <= x <= 3/4 && y == 0]};

Question
--------

> How to compute the gravitational potential created by a disc (and a ring) using FEM in NDSolve?

In a broader sense I think I am asking how can FEM methods deal with PDEs with boundaries at infinity?
*I am guessing that one strategy might be to move the boundary sufficiently far away and increasing sampling within the inner region?*

Note that my attempt above is imposing fixed potential on the disk not fixed density. I am not sure this is important or not, but ideally (to compare to the analytical solution below) fixing density would be better.

**PostScriptum**

I have found this (nice!) [blog] which provides me with an analytic solution
as follows

PhiDiskData =
WolframAlpha[
"electric potential of a charged disk", {{"Result", 1},
"Input"}] // ReleaseHold;
PhiDisk = PhiDiskData /.  QuantityVariable[a_, _] -> a /. {
Q -> Pi R^2, "ElectricConstant" -> 1};
Phi = PhiDisk /. { x -> r Cos[Theta], y -> r Sin[Theta]} //
Simplify[#, Assumptions -> {r > 0}] &

![Mathematica graphics](http://i.stack.imgur.com/N0l1d.png)

Clear[fD]; fD =
FullSimplify[-D[PhiDisk, {{x, y, z}}]  /. x^2 + y^2 -> r^2,
Element[z, Reals] && r > 0]  /. {x -> r Cos[Theta], y -> r Sin[Theta]};
fD = -{Sqrt[fD[]^2 + fD[]^2] // FullSimplify, fD[]};

![Mathematica graphics](http://i.stack.imgur.com/26Cjq.png)

So that

phiN = (Phi /. { Theta -> 0, R -> 1/2}); pl1 =
ContourPlot[Evaluate[phiN], {r, 0, 2}, {z, -2, 2}, Exclusions -> {},
Contours -> 15,ColorFunction -> (ColorData["RedBlueTones"][1 - #] &),
Epilog -> {Thickness[0.02], Line[{{0, 0}, { 1/2, 0}}]},
FrameLabel -> {r, z}, PlotRange -> All, AspectRatio -> 2];
pl2 = StreamPlot[(fD /. R -> 1/2), {r, 0, 2}, {z, -2, 2},
AspectRatio -> 2, StreamStyle -> White];
pl3 = Show[pl1, pl2, PlotRange -> {{0, 1.5}, {-0.5, 1}},
AspectRatio -> 1]

yields

![Mathematica graphics](http://i.stack.imgur.com/bUX03.png)

> So my question amounts to finding this solution numerically.

Note that the analytic solution works nicely for rings as well (if defined as the difference between two discs.)

phiN = (Phi /. { Theta -> 0,
R -> 1}) - (Phi /. { Theta -> 0, R -> 1/2}); pl1 =
ContourPlot[Evaluate[phiN], {r, 0, 2}, {z, -2, 2}, Exclusions -> {},
Contours -> 15,ColorFunction -> (ColorData["RedBlueTones"][1 - #] &),
Epilog -> {Thickness[0.02], Line[{{1/2, 0}, { 1, 0}}]},
FrameLabel -> {r, z}, PlotRange -> All, AspectRatio -> 2];
pl2 = StreamPlot[(fD /. R -> 1) - (fD /. R -> 1/2), {r, 0, 2}, {z, -2,
2},AspectRatio -> 2, StreamStyle -> White];
pl4 = Show[pl1, pl2, PlotRange -> {{0, 1.5}, {-0.5, 1}},
AspectRatio -> 1]

![Mathematica graphics](http://i.stack.imgur.com/RCWcD.png)

: http://blog.wolfram.com/2012/09/27/3d-charges-and-configurations-with-sharp-edges/```