Odd behaviour of the weighting potential

Hi again!

I am unsure the category is correct, but it is what it is I guess.

The problem I am having is that the weighting potential is showing gaps when I plot it. Given that the weighting potential is plotted for pixel (1,1), it seems mostly correct, see the xy-plane figure below.

And looking at the figure in the xz-plane and the yz-plane, the behaviour is mostly correct, except for the parts where the field seems to dissapear at the top, the the two figures below.


We can also see that the weighting field starts >-325um when it should start at -325um exactly, se the figure below, even though it stops at the correct 325um by the electrodes.

So my main question is now

  1. Is this behaviour, both at the top and the bottom, expected due to some odd plotting since it is 0.000 at the white parts, or is this something else?
    I can also add that the fields have been formatted using python and inspected in matplotlib, and all looks good there. And when I output the plots of the E-field from Allpix they also show correctly.

And a second question is
2. Is it hard coded not to plot the full weighting potential, but only for pixel (1,1) and its closest surrounding pixels?

I also have thoughts regarding the extent of the weighting potential. My pixels are highly assymetric (think 10um by few hundred um) as can be seen in the images above.
Since the electrodes are so large in one dimension, and so small in another, I am keen on having the field extending to say 7 pixels in each direction in the small dimension, and only extending 1 in each direction in the large dimension, since the large dimension is much less critical. My question is now

  1. Can I set “distance = 7” in the TransientPropagation module and get the correct behaviour in the small dimension of the electrodes, even though the weighting potential only extends to the neighbouring pixels in the large dimension, or will this cause problems?
    The reason for this is that I want to have a high resolution of the sampling of the fields (or at least the option to have it), but the amount of RAM is limited. And with the charge sharing in the large dimension of the electrodes will be insignificant compared to the charge-sharing in the small dimension, it would be nice to limit the extent of the weighting field there.

Thank you in advance,
Rickard

Hi @rickard,

loads of questions, but all well motivated :smiley: Let me go through them one-by-one:

  1. The white part is just to the value 0.00 not being plotted. This is an artifact from the ROOT plotting I believe.
  2. There seems to be some misconception. There is no such thing as a “full weighting potential” since the weighting potential changes and depends on which electrode you calculate the induced current for. That is why in the calculation of the weighting potential, you always only set one electrode to unit potential and all others to ground. We do something similar in APSQ then: for each electrode (or pixel) we would like to calculate the current for, we center the weighting potential such that the electrode on unit potential lies on the pixel in question. we then calculate the induced current by the carriers moving around somewhere, and then change the position of the WP to calculate the same current for another pixel. See explanation here - additions to make this more understandable would be very welcome!
  3. In the light of answer 2 and looking at your weighting potential, I think “1” should be a good value for you :slight_smile: The weighting potential does not extend even to neighboring pixels to a significant degree, so even in immediatate neighbors you don’t seem to expect significant induced currents. The weighting potential only needs to cover distances from the pixel on unit potential until where you get non-zero values. distance = 7 seems excessive.

Having said that, the weighting potential you are showing loosk suspicious to me, mostly because it has s a curiosly repeating pattern along x but not in y, this does not look like any shape I would expect. Are you sure your boundary conditions are set up properly? Maybe use the WeightingPotential Generator to generate a reference field to compare with - it uses the full pixel size as electride, which seems to be the case also in your plots above.

Best,
Simon

Hi @simonspa,

Thank you for the quick reply!

Regarding 2. I understand that the weighting field is unique for each electrode and that it is translated to each electrode in the simulation, but maybe I formulated it badly. By full weighting field I meant its full spatial width and depth in XY, since it only plots the weighting field for the (1,1) pixel and its extent to one pixel over, whilst the loaded one extends 7 pixels in each directions. I assume that the weighting field is plotted only for the (1,1) pixel and the spatial coordinates in the plot only extends to the neighboring pixels, and that the full extent of the weighting field for the pixel is used in the simulations?

And for 3. The weighting field does not extend significantly, but it is necessary for us to take at least 3-4 in each direction into account. Though I agree that 7 is a bit excessive, but better safe than sorry. And my question remains: in the X-direction, where the pitch is on the order of 10um, compared to Y where the pitch is a few hundred microns, can I have an asymmetric weighting field that extends to 7 neighboring pixels in the X-direction and only the adjacent 1 in the Y-direction?
I assume so since no errors are thrown.
4. In connection to 3. Does APSQ treat the points outside where the weighting field extends to as 0?

And regarding your final comments, could it be due to the highly asymmetric pixels(?) ->The scaling is very off since the pixels are ~504um in Y and 14um in X, with the electrodes being 501umx10um. I have imported the fields into Python and checked, and they seem correct there, centered and all. Though I do use an odd number of points, if that could cause problems(?)
Regarding the solving for the weighting potential, I use COMSOL and solve the Laplace equation with a 15 by 15 pixels grid,

  • all electrodes are at ground (the pixel electrodes and the backside electrode) except the center one,
  • the middle electrode is at unit potential,
  • the sides of the volume have a zero flux condition (Should be fine since it is far far away from the electrode of interest, but input is appreciated).
  • I then export the potential in a regular grid with the electrode of interest in the center.

Do you have any recommendations regarding the process?
I do basically the same for the electric field, but use the physics models of COMSOL.

I really appreciate the help!
/Rickard

Hi @rickard

sorry for the misunderstanding - now things become more clear to me :smiley:

  1. This is correct. The plotting is done for the central pixel and one ring around it, but the full WP will be used in the simulation. It is actually a good idea to extend the plotting to the full WP provided. I opened an issue for this.
  2. This is correct, you can have an asymmetric WP without any issues. What happens in the Y direction is answered in 4. :slightly_smiling_face:
  3. Yes, all points outside the WP provided (along x and y that is) are treated as zero, i.e. there is no contribution in induced current to the pixel in question.

Your procedure for the WP simulation in COMSOL sounds reasonable (having no experience with COMSOL myself). I am not sure what the “zero flux” condition translates to, but there should be a difference between electric field (periodic boundary conditions) and weighting potential (no periodicity). I just wonder why there is a reappearing non-zero contributions at the outside of your outer pixels along x, but that might really be the electrode configuration - and the difference in y could indeed be explained by the order-of-magnitude larger pitch in that coordinate.

Best,
Simon

Hi @simonspa,

No worries at all, I think my explaining skills needs a bit of improvement as well :slight_smile:

Great to hear that I am not completely off with the weighting field, but to expand on it a bit:
The zero flux just says that \overrightarrow{n}\cdot \nabla f = 0 where \overrightarrow{n} is the surface normal, and f is the potential. Does this sound like it would produce the correct weighting potential given the above boundary conditions as well?
The borders of the simulation volume are far away from the electrode of interest, so the effects should be negligible.

But there is something odd with the figures that comes out for the E-field. I don’t have access to the figures now, but I can try to explain: when I plot with output_plots_project = z and output_plots_projection_percentage = 1 I get essentially 0 (just white with nothing plotted) for all figures for all E-components. But if I change it to output_plots_projection_percentage = 0.999 I get some output that looks like what I would expect.
Again, I have double checked in Python that the field is correct, and plotting with output_plots_project = y and output_plots_projection_percentage = 0.5 shows correct. Is this perhaps also some artifact of ROOT plotting?

Best,
Rickard

Hi @rickard

sorry, late response :slightly_smiling_face:

Does this sound like it would produce the correct weighting potential given the above boundary conditions as well?

It does not actually, because for the weighting potential you do expect the potential to cross pixel boundaries. Your zero flux definition however requires that the derivative of the potential is zero, which to me sounds exactly what periodic boundary conditions would do - if you flip the field at the surface, you get a smooth transition with zero flux. But this is not hwo it should be for the WP.

However, you said the borders are far from the electrode, so if your WP is zero there it indeed does not make a difference :slight_smile:

Is this perhaps also some artifact of ROOT plotting?

it might be, there are all sorts of oddities lurking in ROOT :smirk: If you cross-checked your fields I would say: rely on that.

Best,
Simon

Hi @simonspa,

You’re right, I knew it when I posted. I changed to having the surface to floating, and I simulate a pixel grid that is 51 pixels wide to minimize edge effects, and take only the centre 7-11 pixels (overkill I know, but paranoia…), so the edges are far enough away for the field to be essentially 0, so should be ok! I also get good agreement between Allpix and experiments using a sensor prototype!

Another question though, which I thought I would ask here in case anyone else is also wondering. You have implemented the Mesh Converter to accept the DF-ISE generated by Silvaco TCAD. Is there a reason to why you have chosen Silvaco TCAD instead of something else, say COMSOL? Is Silvaco “better” than COMSOL in some ways, or was Silvaco just the software you were used to beforehand?

Best,
Rickard

Hi @rickard

sorry for the delay, we had an issue with the login of the forums for a few days…

Is there a reason to why you have chosen Silvaco TCAD

The sole and only reason we started with Synopsys Sentaurus (and then later added Silvaco) is that these two tools where the ones already in use and with licenses available for us to run.

Is Silvaco “better” than COMSOL in some ways

Not that I am aware of, and there have been a range of users before who have used COMSOL with very nice results.

If you were inclined so, it would be great to provide a direct connector also for COMSOL output files, so we can start closing this gap of support. :slight_smile:

All the best,
Simon

Hi @simonspa,

I see, thank you for the insight!

If I had the knowledge and know-how I would, however, at my current programming level the resulting code would be held together by only hopes and dreams…
For now I will go with translating using Python.

Best,
Rickard

Hi @rickard

feel free to share your Python script with me and I will see to find the time to reimplement it in the mesh converter.

Best,
Simon

Hi @simonspa,

I’m happy to share it. It is not pretty, as I just threw it together and it seems to work.

A peculiar thing I need to mention is that even when making the mesh symmetric in the X- and Y-directions around the center point of the center pixel, the X-component of the E-field at the surface (electrode side) is not symmetric after exporting it and inspecting in Python, even though the mesh is symmetric and the solution is symmetric and proper when inspected in COMSOL. This quickly disapears as one moves into the sensor away from the electrodes.
This is especially odd given that the evaluation points should be in the mirrored points of the mesh nodes.
I am in contact with COMSOL to see if we can figure out why this happens and will try to remember to update this post if a solution comes along.

I really need to stress that it is not a pretty piece of code, and perhaps not suited to be implemented into the framework, especially given the peculiarity above. It would be prettier, if not better, to implement the use of the COMSOL Multiphysics binary file (.mphbin) as defined in the COMSOL programming reference manual, and directly interpolate from the nodes of the solution. Though I suspect that this would need quite a bit more work, by someone much more well versed in C++ than I.

In the mean time, the code is below (I take no responsibility for the code and it is used at any ones own “risk”).

  1. Export the field from COMSOL in a grid, evenly spaced, as a .txt file
  2. Use the code below to format the field into a .init file that APSQ can read
import numpy as np
from tqdm import tqdm

def loadField(path):
    """
    Takes in the file that you want to convert to a format APSQ can read. 
    
    Input
    --------------------------
    path: string variable of the path to the relevant exported text file you get from COMSOL. Can be either the electric field or the weighting potential
    
    Returns a tuple
    --------------------------
    XX, YY, ZZ, field_list,  width, depth, thickness
    

    Input text file formatting
    --------------------------
    The electric field file needs to be formatted as comma separated values, row by row in the format x, y, z, Ex, Ey, Ez.
    The weighting potential file needs to be formatted similarly as x, y, z, V.
    Above the quantities x, y, z is the 3D point where the E-field is evaluated at, and Ex, Ey, Ez is the field vector at that point,
    similarly for the weighting potential, but the weighting potential is a scalar.
    Also, the fields needs to be evaluated at regular intervals, i.e. not in the mesh points of COMSOL.
    
    """
    # with open(path+'weightingField_1um_15wide_3deep_14umPitch.txt', 'r') as f:
    with open(path, 'r') as f:
        lines = f.readlines()
    
    N = len(lines)
    
    n = np.sum([xx[0] == "%" for xx in lines[0:100]]) # In case the header is included
    
    N -= n
    
    x = np.zeros(N)
    y = np.zeros(N)
    z = np.zeros(N)
    
    
    data = lines[100].strip().split(",") #100 is arbitrary, only to get "dim"
    
    # Depending on your formatting in the E field, you might want to also have the potential as the seventh value.
    if (len(data) == 7) or (len(data) == 6):
        dim = 3
    elif len(data) == 4:
        dim = 1
    else:
        raise Exception("Unknown data formatting")

    if dim == 3:
        field3D = np.zeros((N,3))
    else:
        field1D = np.zeros(N)
    
    idx = 0
    lines = tqdm(lines)
    lines.set_description("{: <15}".format("Loading data"))
    for line in lines:
        
        if line[0] == "%":
            continue
        else:
            data = line.strip().split(",")
    
            x[idx] = round(float(data[0]),3)
            y[idx] = round(float(data[1]),3)
            z[idx] = round(float(data[2]),3)
            if dim == 3:
                field3D[idx] = float(data[3]), float(data[4]), float(data[5])
            else:
                field1D[idx] = float(data[3])
            idx += 1
    
    width = round(x.max()-x.min())
    depth = round(y.max()-y.min())
    thickness = round(z.max()-z.min())
    
    diffX = round(np.mean(np.diff(np.unique(x))),4)
    diffY = round(np.mean(np.diff(np.unique(y))),4)
    diffZ = round(np.mean(np.diff(np.unique(z))),4)
    
    #Convert the positions to integers starting from 1, Allpix use them as indices in its arrays
    
    x = (x - x.min())/diffX + 1
    y = (y - y.min())/diffY + 1
    z = (z - z.min())/diffZ + 1
    
    Nx = np.unique(x).size
    Ny = np.unique(y).size
    Nz = np.unique(z).size
    
    if dim == 1:
        field1D[field1D<1e-12] = 0
        field1D[field1D > 1-1e-12] = 1
    
    if dim == 3:
        eFieldX = np.zeros((Ny,Nx,Nz))
        eFieldY = np.zeros((Ny,Nx,Nz))
        eFieldZ = np.zeros((Ny,Nx,Nz))
    else:
        weightingField = np.zeros((Ny,Nx,Nz))
    
    XX = np.array([int(k) for k in x])
    YY = np.array([int(k) for k in y])
    ZZ = np.array([int(k) for k in z])
    
    posRange = tqdm(range(x.size))
    posRange.set_description("{: <15}".format("Filling arrays"))
    for idx in posRange:
        
        xCor, yCor, zCor = XX[idx]-1, YY[idx]-1, ZZ[idx]-1
        if dim == 3:
            eFieldX[yCor, xCor, zCor] = field3D[idx,0]
            eFieldY[yCor, xCor, zCor] = field3D[idx,1]
            eFieldZ[yCor, xCor, zCor] = field3D[idx,2]
        else:
            weightingField[yCor, xCor, zCor] = field1D[idx]
    
    if dim == 3:
        return (XX, YY, ZZ, [eFieldX, eFieldY, eFieldZ],  width, depth, thickness)
    else: 
        return (XX, YY, ZZ, [weightingField],  width, depth, thickness)
    
    
def write_field(file_name, data_tuple):
    """
    This formats the field from COMSOL according to this page: https://allpix-squared.docs.cern.ch/docs/14_additional/mesh_converter/
    """
    (XX, YY, ZZ, field_list, width, depth, thickness) = data_tuple
    dim = len(field_list)
    
    if dim == 3:
        eFieldX, eFieldY, eFieldZ = field_list[0], field_list[1], field_list[2]
    elif dim == 1:
        weightingField = field_list[0]
    
    # Need these as indices for the arrays
    Xs = np.unique(XX)
    Ys = np.unique(YY)
    Zs = np.unique(ZZ)
    
    Nx = Xs.shape[0]
    Ny = Ys.shape[0]
    Nz = Zs.shape[0]
    
    #Check for nans
    if dim == 3:
        nans = np.sum(np.isnan(eFieldX))
        nans += np.sum(np.isnan(eFieldY))
        nans += np.sum(np.isnan(eFieldZ))
    else:
        nans = np.sum(np.isnan(weightingField))
        
    assert nans == 0   
    

    with open(file_name, 'wb') as f:
        X = tqdm(Xs)
        X.set_description("{: <15}".format("Formatting data"))
        line = "Solved_at_250V_293K \n\
    ##SEED##  ##EVENTS## \n\
    ##TURN## ##TILT## 1.0 \n\
    0.00 0.0 0.00 \n\
    {} {} {} 293 0 0 0 {} {} {} 0\n".format(thickness, width, depth, Nx, Ny, Nz).encode('ascii')
        f.write(line)
        for xx in X:
            for yy in Ys:
                for zz in Zs:
                    if dim == 3:
                        line = "{} {} {} {} {} {}\n".format(int(xx), 
                                                            int(yy), 
                                                            int(zz), 
                                                            eFieldX[yy-1,xx-1,zz-1]/100,
                                                            eFieldY[yy-1,xx-1,zz-1]/100,
                                                            eFieldZ[yy-1,xx-1,zz-1]/100).encode('ascii')
                    else:
                        line = "{} {} {} {}\n".format(int(xx), int(yy), int(zz), weightingField[yy-1,xx-1,zz-1]).encode('ascii')
                    f.write(line)



if __name__ == "__main__":
    
    file = 'E_field_path.txt'
    
    # These are just to check if it is a weighting potential or an e-field, change to anything else
    if "E_field" in file:
        data_tuple = loadField(file)
        eFieldX, eFieldY, eFieldZ = data_tuple[3][0], data_tuple[3][1], data_tuple[3][2]
    elif "W_potential" in file:
        data_tuple = loadField(file)
        weightingField = data_tuple[3][0]
    else:
        raise Exception()
    
    dim = len(data_tuple[3])
    
    width, depth, thickness = data_tuple[-3:]
    
    print("\n\nDOUBLE CHECK")
    print("Width: {}\nDepth: {}\nThickness: {}".format(width, depth, thickness))
    
    
    
    if dim == 3:
        fileName = "E_field.init"
    else:
        fileName = "W_potential.init"
    write_field(fileName, data_tuple)
1 Like

Thanks @rickard

do you by any chance have a sample mphbin file for me along with a description of its content?

/Simon

Hi @simonspa,

Sure! I can’t attach them here on the forum, but I will send them to your email.

/Rickard