If you enable RI technique in ORCA, the DFT calculation speed will be significantly faster than Gaussian. It worth to note that Multiwfn is able to easily generate input file of ORCA, namely loading a file into Multiwfn, and then inputting oi, then you can select a calculation level to generate corresponding ORCA input file containing well-suited keywords.

Note that ORCA doesn't support wB97XD but supports its successor, namely wB97X-D3. In ORCA, it is best to use def2-SVP instead of 6-31G**, because there is a built-in corresponding auxiliary basis set in ORCA.

Overall, I recommend you using below ORCA input file for realizing your present purpose.

! wB97X-D3 def2-SVP def2/J RIJCOSX noautostart miniprint nopop
%maxcore  1000
%pal nprocs   4 end
* xyz   0   1
S     -0.09869900    1.71298500   -0.12179500
O     -1.39721400    1.77749800   -0.71589700
O     -0.03593400   -2.04116900    4.37062200
O      0.53563000    2.92178000    0.32404700
N      0.90863500    1.06969100   -1.24145200
C     -0.14139000    0.60412700    1.23827400
C     -1.12393300   -0.35130300    1.33423900
C     -1.10577100   -1.26369100    2.36898600
C     -0.11108300   -1.20009300    3.31568100
C      0.87269900   -0.22660900    3.22619400
C      0.85230200    0.66647900    2.18581100
C      0.51199300   -0.05987300   -2.02410100
C     -0.38281700    0.10303900   -3.05843300
C     -0.71854800   -0.95871000   -3.84651700
C     -0.15123300   -2.19332700   -3.61222500
C      0.74954700   -2.36045400   -2.59060700
C      1.07583800   -1.30258900   -1.78434200
H     -0.66260900   -2.57033700    4.34322000
H      1.68653900    1.41556800   -1.36125300
H     -1.80285500   -0.38407500    0.69872400
H     -1.76370700   -1.91878100    2.42654400
H      1.54236500   -0.17914500    3.87047500
H      1.51372500    1.31786800    2.12141400
H     -0.75823100    0.93831600   -3.21974900
H     -1.32936300   -0.85071100   -4.53977600
H     -0.37966000   -2.91775600   -4.15232100
H      1.13791900   -3.19277500   -2.44374600
H      1.67268500   -1.41911300   -1.08219700
 *

On my Intel 4-core notebook, the calculation finished within 2 minutes.

I have tried using Orca but I find it slow. My question is could Multiwfn be used for this purpose?.

The contents of a typical GJF that runs on G16W in about 5 minutes is below.

%nprocshared=  4
#WB97XD/6-31G** guess=core nosym density=current pop=esp cube=(cards,frozencore)

#hps1xxxx 'P 21 21   fragm   1

0 1
S   -0.098699    1.712985   -0.121795
O   -1.397214    1.777498   -0.715897
O   -0.035934   -2.041169    4.370622
O    0.535630    2.921780    0.324047
N    0.908635    1.069691   -1.241452
C   -0.141390    0.604127    1.238274
C   -1.123933   -0.351303    1.334239
C   -1.105771   -1.263691    2.368986
C   -0.111083   -1.200093    3.315681
C    0.872699   -0.226609    3.226194
C    0.852302    0.666479    2.185811
C    0.511993   -0.059873   -2.024101
C   -0.382817    0.103039   -3.058433
C   -0.718548   -0.958710   -3.846517
C   -0.151233   -2.193327   -3.612225
C    0.749547   -2.360454   -2.590607
C    1.075838   -1.302589   -1.784342
H   -0.662609   -2.570337    4.343220
H    1.686539    1.415568   -1.361253
H   -1.802855   -0.384075    0.698724
H   -1.763707   -1.918781    2.426544
H    1.542365   -0.179145    3.870475
H    1.513725    1.317868    2.121414
H   -0.758231    0.938316   -3.219749
H   -1.329363   -0.850711   -4.539776
H   -0.379660   -2.917756   -4.152321
H    1.137919   -3.192775   -2.443746
H    1.672685   -1.419113   -1.082197

c:\clp\doc\dens\hps1xxxxa.den
    0   -4.760000   -6.279999   -7.080000
  120    0.080000    0.000000    0.000000
  160    0.000000    0.080000    0.000000
  180    0.000000    0.000000    0.080000