from ngsolve import *
from netgen.geom2d import SplineGeometry

#--------------------------------------------------------------------
def TopRectangle(geom):
    
    # Make top subdomain (material 1)  \Omega_1: 
    #          +-------------------------------------------+ (2,2)
    #          |          	                               |
    #          |       Subdomain 1                         |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          +-------------------------------------------+
    #    (-2,-0.35)                                      (2,-0.35)


    # First, make a list of all points:
    #             point 0,       point 1,       point 2,       point 3
    top_pnts = [ (+2.00,-0.35), (+2.00,+2.00), (-2.00,+2.00), (-2.00,-0.35)]
        
    # Add them to given geometry and collect the vertex numbers
    top_nums = [geom.AppendPoint(*p) for p in top_pnts]

    # Make a list of line segments in this format:
    #           beginning point number,
    #           |            end point number,
    #           |            |             boundary condition number,
    #           |            |             |     domain on left side,
    #           |            |             |     |   domain on right side
    #           |            |             |     |   |
    #           V            V             V     V   V
    lines  = [ (top_nums[0], top_nums[1],  10,   1,  0),
               (top_nums[1], top_nums[2],  10,   1,  0),
               (top_nums[2], top_nums[3],  10,   1,  0),
               (top_nums[3], top_nums[0],  10,   1,  0) ]
                
    for p1,p2,bn,dl,dr  in  lines:
        geom.Append([ "line", p1, p2 ],
                    bc=bn, leftdomain=dl, rightdomain=dr )
        
    return geom
#--------------------------------------------------------------------

def PunchCircle(geom):
    #
    #          +-------------------------------------------+ 
    #          |         	                               |
    #          |       Subdomain 1                         |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                                           |
    #          |                3   2   1                  |
    #          |                   ---      Disk of radius |
    #          |                4(  o  )0   0.25 centered  |
    #          |                   ---      at o=origin    |
    #          |                5   6   7                  |
    #          +-------------------------------------------+ 

    disc_pnts = [(+0.25, +0.00),   # point 0
                 (+0.25, +0.25),   # point 1
                 (+0.00, +0.25),   # point 2
                 (-0.25, +0.25),   # point 3
                 (-0.25, +0.00),   # point 4
                 (-0.25, -0.25),   # point 5
                 (+0.00, -0.25),   # point 6
                 (+0.25, -0.25) ]  # point 7

    disc_pnums = [geom.AppendPoint(*p) for p in disc_pnts]

    # Make a SPLINE curve in this tuple format:
    #   (beginning control point, middle control point, final control point,
    #       boundary condition number, domain on left side,  domain on right side)
    #
    curves = [ (disc_pnums[0], disc_pnums[1], disc_pnums[2],  1,0,1),
               (disc_pnums[2], disc_pnums[3], disc_pnums[4],  1,0,1),
               (disc_pnums[4], disc_pnums[5], disc_pnums[6],  1,0,1),
               (disc_pnums[6], disc_pnums[7], disc_pnums[0],  1,0,1)  ]
        
    for p0,p1,p2,bc,left,right in curves:
        geom.Append( ["spline3", p0,p1,p2],
                     bc=bc, leftdomain=left, rightdomain=right)
    return geom
#--------------------------------------------------------------------

geometry = SplineGeometry()
geometry = TopRectangle(geometry)
geometry = PunchCircle(geometry)  
mesh = Mesh( geometry.GenerateMesh() )