Apply a disturbance in IEEE39

Hello everyone,

I want to put a composite load at one bus, for example, bus 20, and apply a 3-phase fault at bus 6 of the IEEE39 system.

I am new to python and PSSE. I have tried to find the data and the model. I have found this file on this website. Since I can not share link, this is the question which I have found the IEEE39 data: "IEEE 39 BUS DYNAMIC SIMULATION IN PSSE" from Jconto's google drive.

I have tried to use the mentioned code and model. However, I have some problems regarding that. 1- I have to change the time of the simulation run, but for all the run time, the outputs channels have a length of 85.

2- I don't have any clue how to add a load to one specific bus (for example a composite load to bus 20).

3- In the code, there is a fault for one line. But how to apply a 3-phase fault?

I have provided the code as follow:

import os, sys
import numpy as np
import csv
import matplotlib.pyplot as plt
sys_path_PSSE = 'C:\Program Files (x86)\PTI\PSSEXplore34\PSSPY27'  # or where else you find the psspy.pyc
sys.path.append(sys_path_PSSE)
os_path_PSSE = 'C:\Program Files (x86)\PTI\PSSEXplore34\PSSBIN'  # or where else you find the psse.exe
os.environ['PATH'] = ';' + os_path_PSSE
import psspy
import dyntools
import redirect
import os, sys
def runscript(myscript):
    import psspy
    if os.path.isfile(myscript):
       root,ext = os.path.splitext(myscript)
   #print root,ext
       if   ext.upper()=='.PY':
            execfile(myscript)
       elif ext.upper()=='.IDV':
            psspy.runrspnsfile(myscript)
       elif ext.upper()=='.IRF':
            psspy.runiplanfile(myscript)
try:
  psspy
  GUI = True
  if sys.argv[-1]=='' or len(sys.argv[-1].strip())<1:
     sys.argv.pop()
except:
  GUI = False
  import psspy
  psspy.psseinit(150000)
_i = psspy.getdefaultint()
_f = psspy.getdefaultreal()
_s = psspy.getdefaultchar()
import redirect
try:
  redirect.psse2py()
except:
  pass
keyword = 'IEEE39'
savkey  = 'IEEE39_v33'
dyrfile = 'IEEE39.dyr'
study   = 'IEEE39_flat'
conl    = 'Conl.idv'
cnvfile = '%s_cnv.sav'%study
channels= 'channels.idv'
snpfile = '%s'%keyword
outfile = '%s.out'%study
logfile = '%s.log'%study
sys.stdout = open(logfile,'w') # redirect all prints to this log file
psspy.case(savkey)
runscript(conl)
psspy.fnsl((_i,0,_i,_i,_i,_i,_i,0))
psspy.cong(0)
psspy.ordr(0)
psspy.fact()
psspy.tysl(0)
psspy.tysl(0)
psspy.save(cnvfile)
psspy.dyre_new([1,1,1,1],dyrfile,"","","")
psspy.snap([-1,-1,-1,-1,-1],snpfile)
sid    = 2       # -1 = All buses, or enter 0-9 for a subsystem:
region = [1]
ierr = psspy.zsys(sid, len(region), region)
print 'psspy.asys - ',ierr
psspy.chsb(2,0,[-1,-1,-1,5,0,0])   # system totals by zone subsystem
psspy.dynamics_solution_param_2([_i, _i, _i, _i, _i, _i, _i, _i], [_f, _f, 0.00001, _f, _f, _f, _f, _f])
psspy.set_relang(1,38,r"""1""")

psspy.strt_2([0,1],outfile)



psspy.run(0,53,990,9,0)
sys.stdout.close()                # ordinary file object
sys.stdout = sys.__stdout__

ch = dyntools.CHNF(outfile)
short_title, chanid_dict, chandata_dict = ch ...