Skip to content
Extraits de code Groupes Projets
semiconductor.html 24,2 ko
Newer Older
  • Learn to ignore specific revisions
  • Nicolas Roisin's avatar
    Nicolas Roisin a validé
    <!doctype html>
    <html lang="en">
    <head>
    <meta charset="utf-8">
    <meta name="viewport" content="width=device-width, initial-scale=1, minimum-scale=1">
    <meta name="generator" content="pdoc3 0.11.6">
    <title>dopes.data_analysis.semiconductor API documentation</title>
    <meta name="description" content="">
    <link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/10up-sanitize.css/13.0.0/sanitize.min.css" integrity="sha512-y1dtMcuvtTMJc1yPgEqF0ZjQbhnc/bFhyvIyVNb9Zk5mIGtqVaAB1Ttl28su8AvFMOY0EwRbAe+HCLqj6W7/KA==" crossorigin>
    <link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/10up-sanitize.css/13.0.0/typography.min.css" integrity="sha512-Y1DYSb995BAfxobCkKepB1BqJJTPrOp3zPL74AWFugHHmmdcvO+C48WLrUOlhGMc0QG7AE3f7gmvvcrmX2fDoA==" crossorigin>
    <link rel="stylesheet" href="https://cdnjs.cloudflare.com/ajax/libs/highlight.js/11.9.0/styles/default.min.css" crossorigin>
    <style>:root{--highlight-color:#fe9}.flex{display:flex !important}body{line-height:1.5em}#content{padding:20px}#sidebar{padding:1.5em;overflow:hidden}#sidebar > *:last-child{margin-bottom:2cm}.http-server-breadcrumbs{font-size:130%;margin:0 0 15px 0}#footer{font-size:.75em;padding:5px 30px;border-top:1px solid #ddd;text-align:right}#footer p{margin:0 0 0 1em;display:inline-block}#footer p:last-child{margin-right:30px}h1,h2,h3,h4,h5{font-weight:300}h1{font-size:2.5em;line-height:1.1em}h2{font-size:1.75em;margin:2em 0 .50em 0}h3{font-size:1.4em;margin:1.6em 0 .7em 0}h4{margin:0;font-size:105%}h1:target,h2:target,h3:target,h4:target,h5:target,h6:target{background:var(--highlight-color);padding:.2em 0}a{color:#058;text-decoration:none;transition:color .2s ease-in-out}a:visited{color:#503}a:hover{color:#b62}.title code{font-weight:bold}h2[id^="header-"]{margin-top:2em}.ident{color:#900;font-weight:bold}pre code{font-size:.8em;line-height:1.4em;padding:1em;display:block}code{background:#f3f3f3;font-family:"DejaVu Sans Mono",monospace;padding:1px 4px;overflow-wrap:break-word}h1 code{background:transparent}pre{border-top:1px solid #ccc;border-bottom:1px solid #ccc;margin:1em 0}#http-server-module-list{display:flex;flex-flow:column}#http-server-module-list div{display:flex}#http-server-module-list dt{min-width:10%}#http-server-module-list p{margin-top:0}.toc ul,#index{list-style-type:none;margin:0;padding:0}#index code{background:transparent}#index h3{border-bottom:1px solid #ddd}#index ul{padding:0}#index h4{margin-top:.6em;font-weight:bold}@media (min-width:200ex){#index .two-column{column-count:2}}@media (min-width:300ex){#index .two-column{column-count:3}}dl{margin-bottom:2em}dl dl:last-child{margin-bottom:4em}dd{margin:0 0 1em 3em}#header-classes + dl > dd{margin-bottom:3em}dd dd{margin-left:2em}dd p{margin:10px 0}.name{background:#eee;font-size:.85em;padding:5px 10px;display:inline-block;min-width:40%}.name:hover{background:#e0e0e0}dt:target .name{background:var(--highlight-color)}.name > span:first-child{white-space:nowrap}.name.class > span:nth-child(2){margin-left:.4em}.inherited{color:#999;border-left:5px solid #eee;padding-left:1em}.inheritance em{font-style:normal;font-weight:bold}.desc h2{font-weight:400;font-size:1.25em}.desc h3{font-size:1em}.desc dt code{background:inherit}.source > summary,.git-link-div{color:#666;text-align:right;font-weight:400;font-size:.8em;text-transform:uppercase}.source summary > *{white-space:nowrap;cursor:pointer}.git-link{color:inherit;margin-left:1em}.source pre{max-height:500px;overflow:auto;margin:0}.source pre code{font-size:12px;overflow:visible;min-width:max-content}.hlist{list-style:none}.hlist li{display:inline}.hlist li:after{content:',\2002'}.hlist li:last-child:after{content:none}.hlist .hlist{display:inline;padding-left:1em}img{max-width:100%}td{padding:0 .5em}.admonition{padding:.1em 1em;margin:1em 0}.admonition-title{font-weight:bold}.admonition.note,.admonition.info,.admonition.important{background:#aef}.admonition.todo,.admonition.versionadded,.admonition.tip,.admonition.hint{background:#dfd}.admonition.warning,.admonition.versionchanged,.admonition.deprecated{background:#fd4}.admonition.error,.admonition.danger,.admonition.caution{background:lightpink}</style>
    <style media="screen and (min-width: 700px)">@media screen and (min-width:700px){#sidebar{width:30%;height:100vh;overflow:auto;position:sticky;top:0}#content{width:70%;max-width:100ch;padding:3em 4em;border-left:1px solid #ddd}pre code{font-size:1em}.name{font-size:1em}main{display:flex;flex-direction:row-reverse;justify-content:flex-end}.toc ul ul,#index ul ul{padding-left:1em}.toc > ul > li{margin-top:.5em}}</style>
    <style media="print">@media print{#sidebar h1{page-break-before:always}.source{display:none}}@media print{*{background:transparent !important;color:#000 !important;box-shadow:none !important;text-shadow:none !important}a[href]:after{content:" (" attr(href) ")";font-size:90%}a[href][title]:after{content:none}abbr[title]:after{content:" (" attr(title) ")"}.ir a:after,a[href^="javascript:"]:after,a[href^="#"]:after{content:""}pre,blockquote{border:1px solid #999;page-break-inside:avoid}thead{display:table-header-group}tr,img{page-break-inside:avoid}img{max-width:100% !important}@page{margin:0.5cm}p,h2,h3{orphans:3;widows:3}h1,h2,h3,h4,h5,h6{page-break-after:avoid}}</style>
    <script defer src="https://cdnjs.cloudflare.com/ajax/libs/highlight.js/11.9.0/highlight.min.js" integrity="sha512-D9gUyxqja7hBtkWpPWGt9wfbfaMGVt9gnyCvYa+jojwwPHLCzUm5i8rpk7vD7wNee9bA35eYIjobYPaQuKS1MQ==" crossorigin></script>
    <script>window.addEventListener('DOMContentLoaded', () => {
    hljs.configure({languages: ['bash', 'css', 'diff', 'graphql', 'ini', 'javascript', 'json', 'plaintext', 'python', 'python-repl', 'rust', 'shell', 'sql', 'typescript', 'xml', 'yaml']});
    hljs.highlightAll();
    /* Collapse source docstrings */
    setTimeout(() => {
    [...document.querySelectorAll('.hljs.language-python > .hljs-string')]
    .filter(el => el.innerHTML.length > 200 && ['"""', "'''"].includes(el.innerHTML.substring(0, 3)))
    .forEach(el => {
    let d = document.createElement('details');
    d.classList.add('hljs-string');
    d.innerHTML = '<summary>"""</summary>' + el.innerHTML.substring(3);
    el.replaceWith(d);
    });
    }, 100);
    })</script>
    </head>
    <body>
    <main>
    <article id="content">
    <header>
    <h1 class="title">Module <code>dopes.data_analysis.semiconductor</code></h1>
    </header>
    <section id="section-intro">
    </section>
    <section>
    </section>
    <section>
    </section>
    <section>
    <h2 class="section-title" id="header-functions">Functions</h2>
    <dl>
    <dt id="dopes.data_analysis.semiconductor.intrinsic_concentration"><code class="name flex">
    <span>def <span class="ident">intrinsic_concentration</span></span>(<span>temp)</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def intrinsic_concentration(temp):
        &#34;&#34;&#34; Function to calculate the intrinsic concentration of silicon from K. Misiakos and Tsamakis, D., “Accurate measurements of the silicon intrinsic carrier density from 78 to 340 K”, Journal of Applied Physics, vol. 74, no. 5, p. 3293, 1993.
        
            args:
                - temp (scalar): the temperature 
                
            return:
                - ni (scalar): the intrinsic concentration in silicon
        &#34;&#34;&#34;  
        
        return 5.29e19 * (temp/300)**2.54 * np.exp(-6726/temp)</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the intrinsic concentration of silicon from K. Misiakos and Tsamakis, D., “Accurate measurements of the silicon intrinsic carrier density from 78 to 340 K”, Journal of Applied Physics, vol. 74, no. 5, p. 3293, 1993.</p>
    <p>args:
    - temp (scalar): the temperature </p>
    <p>return:
    - ni (scalar): the intrinsic concentration in silicon</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.mobility_impurity"><code class="name flex">
    <span>def <span class="ident">mobility_impurity</span></span>(<span>mu_0, carrier='n', Ni=1000000000000000.0, temp=300, dopant='phosphorus')</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def mobility_impurity(mu_0,carrier=&#34;n&#34;,Ni=1e15,temp=300, dopant=&#34;phosphorus&#34;):
        &#34;&#34;&#34; Function to calculate the silicon mobility according to Masetti relation (1983)
        
            args:
                - carrier (string): &#34;n&#34; for electrons, &#34;p&#34; for holes
                - temp (scalar): the temperature 
                - Ni (scalar): the impurity cencentration in cm-3
                - dopant (string): the type of n-type dopant. &#34;phosphorus&#34; or &#34;arsenic&#34;
                
            return:
                - mu_LI (scalar): the electron or hole mobility with the impurity scattering taken into account
        &#34;&#34;&#34;  
        
        # Values are taken from Masetti et al. (1983)
        if dopant==&#34;phosphorus&#34;:
            param_300_n={&#34;mu_min&#34;:68.5,&#34;Cref&#34;:9.20e16,&#34;alpha&#34;:0.711} 
            correction_n={&#34;mu_min&#34;:56.1,&#34;Cref&#34;:3.41e20,&#34;alpha&#34;:1.98} 
            # mu_0=1414
        elif dopant==&#34;arsenic&#34;:
            param_300_n={&#34;mu_min&#34;:52.2,&#34;Cref&#34;:9.68e16,&#34;alpha&#34;:0.680} 
            correction_n={&#34;mu_min&#34;:43.4,&#34;Cref&#34;:3.43e20,&#34;alpha&#34;:2.00} 
            # mu_0=1417
    
        param_300_p={&#34;mu_min&#34;:44.9,&#34;Cref&#34;:22.3e16,&#34;alpha&#34;:0.72}
        correction_p={&#34;mu_min&#34;:29.0,&#34;Cref&#34;:6.1e20,&#34;alpha&#34;:2.0}
        
        expon_temp={&#34;mu_min&#34;:-0.45,&#34;Cref&#34;:3.2,&#34;alpha&#34;:0.065}
         
        if carrier==&#34;n&#34;:
            param_300=param_300_n
            correction=correction_n[&#34;mu_min&#34;]/(1+(correction_n[&#34;Cref&#34;]/Ni)**correction_n[&#34;alpha&#34;])
        else:
            param_300=param_300_p
            correction=correction_p[&#34;mu_min&#34;]/(1+(correction_p[&#34;Cref&#34;]/Ni)**correction_p[&#34;alpha&#34;])
            # mu_0=470.5
    
    
        mu_min=param_300[&#34;mu_min&#34;]*(temp/300)**expon_temp[&#34;mu_min&#34;]
        Cref=param_300[&#34;Cref&#34;]*(temp/300)**expon_temp[&#34;Cref&#34;]
        alpha=param_300[&#34;alpha&#34;]*(temp/300)**expon_temp[&#34;alpha&#34;]
        
        mu_LI=mu_min+(mu_0-mu_min) / ( 1 + ( Ni / Cref )**alpha )-correction
        
        return mu_LI</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the silicon mobility according to Masetti relation (1983)</p>
    <p>args:
    - carrier (string): "n" for electrons, "p" for holes
    - temp (scalar): the temperature
    - Ni (scalar): the impurity cencentration in cm-3
    - dopant (string): the type of n-type dopant. "phosphorus" or "arsenic"</p>
    <p>return:
    - mu_LI (scalar): the electron or hole mobility with the impurity scattering taken into account</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.piezo_ratio_temperature"><code class="name flex">
    <span>def <span class="ident">piezo_ratio_temperature</span></span>(<span>temp, carrier='n')</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def piezo_ratio_temperature(temp,carrier=&#34;n&#34;):
        &#34;&#34;&#34; Function to calculate the correcting ratio due to temperature variations for the piezoresistive coefficients. The data are taken from Kanda (1982) between 200 K and 400K.
            The correction is refered to the 300K coefficients. The new piezoresistive coefficients can be calculated by piezo_coefficient(temp) = ratio * piezo_coefficient(300K)
            
            args:
                - temp (scalar or sequence): the temperatures for which the correting ratio has to be calculated
                - carrier (string): &#34;n&#34; for electrons, &#34;p&#34; for holes
                    
            return:
                - a scalar or sequence with the same dimension as temp with the correcting ratio
        &#34;&#34;&#34;  
    
        temp_kanda1982=np.linspace(125, -75,9)+273
        dpi_kanda1982_n=np.array([0.7547318611987381,0.8067823343848579,0.8611987381703468,0.9298107255520504,1.0007886435331228,1.0977917981072554,1.2113564668769716,1.3438485804416402,1.5165615141955835])
        dpi_kanda1982_p=np.array([0.7523342983957861,0.8037360206464299,0.8598127702627387,0.922894558591104,0.995324258126102,1.0934548187864221,1.205601887444947,1.3411219494071012,1.5093457676819353])
        
        fn=interp.interp1d(temp_kanda1982,dpi_kanda1982_n,kind=&#34;linear&#34;,fill_value=&#34;extrapolate&#34; )
        fp=interp.interp1d(temp_kanda1982,dpi_kanda1982_p,kind=&#34;linear&#34;,fill_value=&#34;extrapolate&#34; )
    
        if carrier==&#34;n&#34;:
            return fn(temp)
        elif carrier ==&#34;p&#34;:
            return fp(temp)
        else:
            return 0</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the correcting ratio due to temperature variations for the piezoresistive coefficients. The data are taken from Kanda (1982) between 200 K and 400K.
    The correction is refered to the 300K coefficients. The new piezoresistive coefficients can be calculated by piezo_coefficient(temp) = ratio * piezo_coefficient(300K)</p>
    <p>args:
    - temp (scalar or sequence): the temperatures for which the correting ratio has to be calculated
    - carrier (string): "n" for electrons, "p" for holes</p>
    <p>return:
    - a scalar or sequence with the same dimension as temp with the correcting ratio</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.piezoresistivity_strain"><code class="name flex">
    <span>def <span class="ident">piezoresistivity_strain</span></span>(<span>strain_tensor, pi11, pi12, pi44)</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def piezoresistivity_strain(strain_tensor, pi11, pi12, pi44):
        &#34;&#34;&#34; Function to calculate the relative change due to stress in silicon
            
            args:
                - stress_tensor (numpy array): the stress tensor for which the stress should be calculated. The voigt notation should be used with a 1 x 6 vector but the function can handle  3 x 3 matrix but only take the upper half in this case.
                - pi11, pi12 and pi44 (scalar): piezoresistive coefficients used to calculate the variations in the relative resistivity. 
                  Values from Smith (1954) are pi11 = -1022 TPa-1, pi12 = 534 TPa-1 and pi44 = -136 TPa-1 for the electrons,
                  and pi11 = 66 TPa-1, pi12 = -11 TPa-1 and pi44 = 1381 TPa-1 for the holes.
                    
            return:
                - an 1 x 6 tensor using the Voigt notation with the relative variation of the resistivity.
        &#34;&#34;&#34;  
    
        piezo_tensor=np.array([[pi11,pi12,pi12,0,0,0],
                                    [pi12,pi11,pi12,0,0,0],
                                    [pi12,pi12,pi11,0,0,0],
                                    [0,0,0,pi44,0,0],
                                    [0,0,0,0,pi44,0],
                                    [0,0,0,0,0,pi44]])
        
        stress_tensor= mec.stress_from_strain(strain_tensor)
        
        return piezo_tensor @ stress_tensor</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the relative change due to stress in silicon</p>
    <p>args:
    - stress_tensor (numpy array): the stress tensor for which the stress should be calculated. The voigt notation should be used with a 1 x 6 vector but the function can handle
    3 x 3 matrix but only take the upper half in this case.
    - pi11, pi12 and pi44 (scalar): piezoresistive coefficients used to calculate the variations in the relative resistivity.
    Values from Smith (1954) are pi11 = -1022 TPa-1, pi12 = 534 TPa-1 and pi44 = -136 TPa-1 for the electrons,
    and pi11 = 66 TPa-1, pi12 = -11 TPa-1 and pi44 = 1381 TPa-1 for the holes.</p>
    <p>return:
    - an 1 x 6 tensor using the Voigt notation with the relative variation of the resistivity.</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.piezoresistivity_stress"><code class="name flex">
    <span>def <span class="ident">piezoresistivity_stress</span></span>(<span>stress_tensor, pi11, pi12, pi44)</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def piezoresistivity_stress(stress_tensor, pi11, pi12, pi44):
        &#34;&#34;&#34; Function to calculate the relative change due to stress in silicon
            
            args:
                - stress_tensor (numpy array): the stress tensor for which the stress should be calculated. The voigt notation should be used with a 1 x 6 vector but the function can handle  3 x 3 matrix but only take the upper half in this case.
                - pi11, pi12 and pi44 (scalar): piezoresistive coefficients used to calculate the variations in the relative resistivity. 
                  Values from Smith (1954) are pi11 = -1022 TPa-1, pi12 = 534 TPa-1 and pi44 = -136 TPa-1 for the electrons,
                  and pi11 = 66 TPa-1, pi12 = -11 TPa-1 and pi44 = 1381 TPa-1 for the holes.
                    
            return:
                - an 1 x 6 tensor using the Voigt notation with the relative variation of the resistivity.
        &#34;&#34;&#34;  
        
        stress_voigt=np.zeros((6,1))
        stress_shape=np.shape(stress_tensor)
        
        if len(stress_shape)==2:
            if stress_shape[0]==3 and stress_shape[1]==3:
                stress_voigt[0]=stress_tensor[0,0]
                stress_voigt[1]=stress_tensor[1,1]
                stress_voigt[2]=stress_tensor[2,2]
                stress_voigt[3]=stress_tensor[1,2]
                stress_voigt[4]=stress_tensor[0,2]
                stress_voigt[5]=stress_tensor[0,1]
            if stress_shape[0]==1 and stress_shape[1]==6:
                stress_voigt=stress_tensor
        else:
            stress_voigt=np.array([stress_tensor])
        
        piezo_tensor=np.array([[pi11,pi12,pi12,0,0,0],
                                    [pi12,pi11,pi12,0,0,0],
                                    [pi12,pi12,pi11,0,0,0],
                                    [0,0,0,pi44,0,0],
                                    [0,0,0,0,pi44,0],
                                    [0,0,0,0,0,pi44]])
        
        return piezo_tensor @ stress_voigt</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the relative change due to stress in silicon</p>
    <p>args:
    - stress_tensor (numpy array): the stress tensor for which the stress should be calculated. The voigt notation should be used with a 1 x 6 vector but the function can handle
    3 x 3 matrix but only take the upper half in this case.
    - pi11, pi12 and pi44 (scalar): piezoresistive coefficients used to calculate the variations in the relative resistivity.
    Values from Smith (1954) are pi11 = -1022 TPa-1, pi12 = 534 TPa-1 and pi44 = -136 TPa-1 for the electrons,
    and pi11 = 66 TPa-1, pi12 = -11 TPa-1 and pi44 = 1381 TPa-1 for the holes.</p>
    <p>return:
    - an 1 x 6 tensor using the Voigt notation with the relative variation of the resistivity.</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.tau_srh"><code class="name flex">
    <span>def <span class="ident">tau_srh</span></span>(<span>Ni, tau_0=5e-05, Nc=5e+16, A=1, B=1, C=0, E=0)</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def tau_srh(Ni, tau_0=5e-5,  Nc=5e16, A=1, B=1, C=0, E=0):
        &#34;&#34;&#34; Function to calculate the SRH lifetime of silicon. The default parameters, similar for electron and hole, are those from D&#39;Avanzo, D. C., Vanzi, M., Dutton, R. W.: One-Dimensional Semiconductor Device Analysis
            (SEDAN). Report G-201-5, Stanford University, 1979.
            
            args:
                - Ni (scalar): the impurity density in the material 
                - tau_0 (scalar): intial value for the lifetime 
                - Nc (scalar): the critical impurity level 
                - A, B, C and D (scalar): the coefficient for the model 
                
            return:
                - the SRH lifetime in silicon
        &#34;&#34;&#34;  
        
        return tau_0/(A+B*(Ni/Nc)+C*(Ni/Nc)**E)</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the SRH lifetime of silicon. The default parameters, similar for electron and hole, are those from D'Avanzo, D. C., Vanzi, M., Dutton, R. W.: One-Dimensional Semiconductor Device Analysis
    (SEDAN). Report G-201-5, Stanford University, 1979.</p>
    <p>args:
    - Ni (scalar): the impurity density in the material
    - tau_0 (scalar): intial value for the lifetime
    - Nc (scalar): the critical impurity level
    - A, B, C and D (scalar): the coefficient for the model </p>
    <p>return:
    - the SRH lifetime in silicon</p></div>
    </dd>
    <dt id="dopes.data_analysis.semiconductor.tau_trap"><code class="name flex">
    <span>def <span class="ident">tau_trap</span></span>(<span>tau_n, tau_p, NA, ND, Etrap=0.56, Eg=1.12, temp=300)</span>
    </code></dt>
    <dd>
    <details class="source">
    <summary>
    <span>Expand source code</span>
    </summary>
    <pre><code class="python">def tau_trap(tau_n, tau_p, NA, ND, Etrap=0.56,Eg=1.12, temp=300):
        &#34;&#34;&#34; Function to calculate the lifetime of silicon due to a trap.
            
            args:
                - tau_n and tau_p (scalar): the life time of the electron and hole due to the trap density. 
                  The lifetime can be calculated by tau_n = 1 / (sigma_n * v_th * Nt ) where sigma_n is the capture cross section, vth is the thermal velocity and Nt is the trap density. 
                  Typical values for the thermal velocities are 2.3e7 and 1.65e7 cm/s for the electrons and holes, respectively. 
                  Typical value for the capture cross section is 1e-15 cm² for a neutral defect. 
                  Typical value for the trap density is 1e12 cm-3 for a neutral defect.  
                - NA and ND (scalar): the acceptor and donor density of the pn junction 
                - Etrap (scalar): the trap level referred to the maximum of the valence band 
                - Eg (scalar): the bandgap of the material 
                - temp (scalar): the temperature 
                
            return:
                - the trap lifetime 
        &#34;&#34;&#34;  
        kB = 1.38e-23 # J/K
        q = 1.602e-19 # C
        
        nieff=intrinsic_concentration(temp)
        
        n_0 = nieff**2 / NA
        p_0 = nieff**2 / ND
        n1 = nieff * np.exp(- q * (Eg - Etrap) / (kB*temp))
        p1 = nieff * np.exp(- q * Etrap / (kB*temp))
        
        return (tau_p * (n_0 + n1) + tau_n * (p_0 + p1)) / (p_0 + n_0)</code></pre>
    </details>
    <div class="desc"><p>Function to calculate the lifetime of silicon due to a trap.</p>
    <p>args:
    - tau_n and tau_p (scalar): the life time of the electron and hole due to the trap density.
    The lifetime can be calculated by tau_n = 1 / (sigma_n * v_th * Nt ) where sigma_n is the capture cross section, vth is the thermal velocity and Nt is the trap density.
    Typical values for the thermal velocities are 2.3e7 and 1.65e7 cm/s for the electrons and holes, respectively.
    Typical value for the capture cross section is 1e-15 cm² for a neutral defect.
    Typical value for the trap density is 1e12 cm-3 for a neutral defect.<br>
    - NA and ND (scalar): the acceptor and donor density of the pn junction
    - Etrap (scalar): the trap level referred to the maximum of the valence band
    - Eg (scalar): the bandgap of the material
    - temp (scalar): the temperature </p>
    <p>return:
    - the trap lifetime</p></div>
    </dd>
    </dl>
    </section>
    <section>
    </section>
    </article>
    <nav id="sidebar">
    <div class="toc">
    <ul></ul>
    </div>
    <ul id="index">
    <li><h3>Super-module</h3>
    <ul>
    <li><code><a title="dopes.data_analysis" href="index.html">dopes.data_analysis</a></code></li>
    </ul>
    </li>
    <li><h3><a href="#header-functions">Functions</a></h3>
    <ul class="">
    <li><code><a title="dopes.data_analysis.semiconductor.intrinsic_concentration" href="#dopes.data_analysis.semiconductor.intrinsic_concentration">intrinsic_concentration</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.mobility_impurity" href="#dopes.data_analysis.semiconductor.mobility_impurity">mobility_impurity</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.piezo_ratio_temperature" href="#dopes.data_analysis.semiconductor.piezo_ratio_temperature">piezo_ratio_temperature</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.piezoresistivity_strain" href="#dopes.data_analysis.semiconductor.piezoresistivity_strain">piezoresistivity_strain</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.piezoresistivity_stress" href="#dopes.data_analysis.semiconductor.piezoresistivity_stress">piezoresistivity_stress</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.tau_srh" href="#dopes.data_analysis.semiconductor.tau_srh">tau_srh</a></code></li>
    <li><code><a title="dopes.data_analysis.semiconductor.tau_trap" href="#dopes.data_analysis.semiconductor.tau_trap">tau_trap</a></code></li>
    </ul>
    </li>
    </ul>
    </nav>
    </main>
    <footer id="footer">
    <p>Generated by <a href="https://pdoc3.github.io/pdoc" title="pdoc: Python API documentation generator"><cite>pdoc</cite> 0.11.6</a>.</p>
    </footer>
    </body>
    </html>